ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

linalg: update to netlib lapack-3.7.1

Browse Source
This commit is contained in:
Christoph Junghans
2018-05-18 15:17:13 -06:00
parent 858c211fdc
commit 2e7b919774
171 changed files with 5320 additions and 4660 deletions
Download Patch File Download Diff File Expand all files Collapse all files

96
lib/linalg/dlasq2.f
View File

@ -2,38 +2,38 @@
*
* =========== DOCUMENTATION ===========
*
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
*> Download DLASQ2 + dependencies
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlasq2.f">
*> [TGZ]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlasq2.f">
*> [ZIP]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlasq2.f">
*> Download DLASQ2 + dependencies
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlasq2.f">
*> [TGZ]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlasq2.f">
*> [ZIP]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlasq2.f">
*> [TXT]</a>
*> \endhtmlonly
*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE DLASQ2( N, Z, INFO )
*
*
* .. Scalar Arguments ..
* INTEGER INFO, N
* ..
* .. Array Arguments ..
* DOUBLE PRECISION Z( * )
* ..
*
*
*
*> \par Purpose:
* =============
*>
*> \verbatim
*>
*> DLASQ2 computes all the eigenvalues of the symmetric positive
*> DLASQ2 computes all the eigenvalues of the symmetric positive
*> definite tridiagonal matrix associated with the qd array Z to high
*> relative accuracy are computed to high relative accuracy, in the
*> absence of denormalization, underflow and overflow.
@ -83,19 +83,19 @@
*> = 2, current block of Z not diagonalized after 100*N
*> iterations (in inner while loop). On exit Z holds
*> a qd array with the same eigenvalues as the given Z.
*> = 3, termination criterion of outer while loop not met
*> = 3, termination criterion of outer while loop not met
*> (program created more than N unreduced blocks)
*> \endverbatim
*
* Authors:
* ========
*
*> \author Univ. of Tennessee
*> \author Univ. of California Berkeley
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*> \author Univ. of Tennessee
*> \author Univ. of California Berkeley
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*
*> \date September 2012
*> \date December 2016
*
*> \ingroup auxOTHERcomputational
*
@ -112,10 +112,10 @@
* =====================================================================
SUBROUTINE DLASQ2( N, Z, INFO )
*
* -- LAPACK computational routine (version 3.4.2) --
* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* September 2012
* December 2016
*
* .. Scalar Arguments ..
INTEGER INFO, N
@ -136,7 +136,7 @@
* .. Local Scalars ..
LOGICAL IEEE
INTEGER I0, I1, I4, IINFO, IPN4, ITER, IWHILA, IWHILB,
$ K, KMIN, N0, N1, NBIG, NDIV, NFAIL, PP, SPLT,
$ K, KMIN, N0, N1, NBIG, NDIV, NFAIL, PP, SPLT,
$ TTYPE
DOUBLE PRECISION D, DEE, DEEMIN, DESIG, DMIN, DMIN1, DMIN2, DN,
$ DN1, DN2, E, EMAX, EMIN, EPS, G, OLDEMN, QMAX,
@ -155,7 +155,7 @@
INTRINSIC ABS, DBLE, MAX, MIN, SQRT
* ..
* .. Executable Statements ..
*
*
* Test the input arguments.
* (in case DLASQ2 is not called by DLASQ1)
*
@ -195,7 +195,7 @@
END IF
Z( 5 ) = Z( 1 ) + Z( 2 ) + Z( 3 )
IF( Z( 2 ).GT.Z( 3 )*TOL2 ) THEN
T = HALF*( ( Z( 1 )-Z( 3 ) )+Z( 2 ) )
T = HALF*( ( Z( 1 )-Z( 3 ) )+Z( 2 ) )
S = Z( 3 )*( Z( 2 ) / T )
IF( S.LE.T ) THEN
S = Z( 3 )*( Z( 2 ) / ( T*( ONE+SQRT( ONE+S / T ) ) ) )
@ -264,19 +264,19 @@
Z( 2*N-1 ) = ZERO
RETURN
END IF
*
*
* Check whether the machine is IEEE conformable.
*
*
IEEE = ILAENV( 10, 'DLASQ2', 'N', 1, 2, 3, 4 ).EQ.1 .AND.
$ ILAENV( 11, 'DLASQ2', 'N', 1, 2, 3, 4 ).EQ.1
*
$ ILAENV( 11, 'DLASQ2', 'N', 1, 2, 3, 4 ).EQ.1
*
* Rearrange data for locality: Z=(q1,qq1,e1,ee1,q2,qq2,e2,ee2,...).
*
DO 30 K = 2*N, 2, -2
Z( 2*K ) = ZERO
Z( 2*K-1 ) = Z( K )
Z( 2*K-2 ) = ZERO
Z( 2*K-3 ) = Z( K-1 )
Z( 2*K ) = ZERO
Z( 2*K-1 ) = Z( K )
Z( 2*K-2 ) = ZERO
Z( 2*K-3 ) = Z( K-1 )
30 CONTINUE
*
I0 = 1
@ -333,7 +333,7 @@
D = Z( I4+1 )*( D / Z( I4-2*PP-2 ) )
END IF
EMIN = MIN( EMIN, Z( I4-2*PP ) )
60 CONTINUE
60 CONTINUE
Z( 4*N0-PP-2 ) = D
*
* Now find qmax.
@ -364,14 +364,14 @@
NDIV = 2*( N0-I0 )
*
DO 160 IWHILA = 1, N + 1
IF( N0.LT.1 )
IF( N0.LT.1 )
$ GO TO 170
*
* While array unfinished do
* While array unfinished do
*
* E(N0) holds the value of SIGMA when submatrix in I0:N0
* splits from the rest of the array, but is negated.
*
*
DESIG = ZERO
IF( N0.EQ.N ) THEN
SIGMA = ZERO
@ -386,7 +386,7 @@
* Find last unreduced submatrix's top index I0, find QMAX and
* EMIN. Find Gershgorin-type bound if Q's much greater than E's.
*
EMAX = ZERO
EMAX = ZERO
IF( N0.GT.I0 ) THEN
EMIN = ABS( Z( 4*N0-5 ) )
ELSE
@ -404,7 +404,7 @@
QMAX = MAX( QMAX, Z( I4-7 )+Z( I4-5 ) )
EMIN = MIN( EMIN, Z( I4-5 ) )
90 CONTINUE
I4 = 4
I4 = 4
*
100 CONTINUE
I0 = I4 / 4
@ -421,7 +421,7 @@
KMIN = ( I4+3 )/4
END IF
110 CONTINUE
IF( (KMIN-I0)*2.LT.N0-KMIN .AND.
IF( (KMIN-I0)*2.LT.N0-KMIN .AND.
$ DEEMIN.LE.HALF*Z(4*N0-3) ) THEN
IPN4 = 4*( I0+N0 )
PP = 2
@ -446,15 +446,15 @@
*
DMIN = -MAX( ZERO, QMIN-TWO*SQRT( QMIN )*SQRT( EMAX ) )
*
* Now I0:N0 is unreduced.
* Now I0:N0 is unreduced.
* PP = 0 for ping, PP = 1 for pong.
* PP = 2 indicates that flipping was applied to the Z array and
* and that the tests for deflation upon entry in DLASQ3
* and that the tests for deflation upon entry in DLASQ3
* should not be performed.
*
NBIG = 100*( N0-I0+1 )
DO 140 IWHILB = 1, NBIG
IF( I0.GT.N0 )
IF( I0.GT.N0 )
$ GO TO 150
*
* While submatrix unfinished take a good dqds step.
@ -497,8 +497,8 @@
140 CONTINUE
*
INFO = 2
*
* Maximum number of iterations exceeded, restore the shift
*
* Maximum number of iterations exceeded, restore the shift
* SIGMA and place the new d's and e's in a qd array.
* This might need to be done for several blocks
*
@ -549,16 +549,16 @@
INFO = 3
RETURN
*
* end IWHILA
* end IWHILA
*
170 CONTINUE
*
*
* Move q's to the front.
*
*
DO 180 K = 2, N
Z( K ) = Z( 4*K-3 )
180 CONTINUE
*
*
* Sort and compute sum of eigenvalues.
*
CALL DLASRT( 'D', N, Z, IINFO )
@ -570,7 +570,7 @@
*
* Store trace, sum(eigenvalues) and information on performance.
*
Z( 2*N+1 ) = TRACE
Z( 2*N+1 ) = TRACE
Z( 2*N+2 ) = E
Z( 2*N+3 ) = DBLE( ITER )
Z( 2*N+4 ) = DBLE( NDIV ) / DBLE( N**2 )