scnc
/
quantum-espresso
mirror of https://gitlab.com/QEF/q-e.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Fixed an array out of bound in Gram-Schmidt. Cleanup of parallel Davidson. 

C.S.


git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@3061 c92efa57-630b-4861-b058-cf58834340f0
This commit is contained in:
sbraccia 2006-05-01 22:14:22 +00:00
parent 98acbc281a
commit 23c45d3ebc
4 changed files with 112 additions and 141 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

117
PW/cegterg.f90
View File

@ -12,7 +12,7 @@
!
!----------------------------------------------------------------------------
SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
overlap, e, btype, notcnv, lrot, dav_iter )
uspp, e, btype, notcnv, lrot, dav_iter )
!----------------------------------------------------------------------------
!
! ... iterative solution of the eigenvalue problem:
@ -20,7 +20,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
! ... ( H - e S ) * evc = 0
!
! ... where H is an hermitean operator, e is a real scalar,
! ... S is an overlap matrix, evc is a complex vector
! ... S is an uspp matrix, evc is a complex vector
!
USE io_global, ONLY : stdout
USE kinds, ONLY : DP
@ -45,7 +45,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
! energy threshold for convergence :
! root improvement is stopped, when two consecutive estimates of the root
! differ by less than ethr.
LOGICAL, INTENT(IN) :: overlap
LOGICAL, INTENT(IN) :: uspp
! if .FALSE. : do not calculate S|psi>
INTEGER, INTENT(IN) :: btype(nvec)
! band type ( 1 = occupied, 0 = empty )
@ -79,8 +79,6 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
! work space, contains psi
! the product of H and psi
! the product of S and psi
COMPLEX(DP) :: eau
! auxiliary complex variable
REAL(DP), ALLOCATABLE :: ew(:)
! eigenvalues of the reduced hamiltonian
LOGICAL, ALLOCATABLE :: conv(:)
@ -161,8 +159,8 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
conv = .FALSE.
!
spsi = ZERO
psi = ZERO
hpsi = ZERO
psi = ZERO
psi(:,:,1:nvec) = evc(:,:,1:nvec)
!
! ... hpsi contains H times the basis vectors
@ -173,7 +171,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
IF ( lelfield ) CALL h_epsi_her( ndmx, ndim, nvec, psi, hpsi )
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL s_psi_nc( ndmx, ndim, nvec, psi, spsi )
!
@ -191,7 +189,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
IF ( lelfield ) CALL h_epsi_her( ndmx, ndim, nvec, psi, hpsi )
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL s_psi( ndmx, ndim, nvec, psi, spsi )
!
@ -214,13 +212,14 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
CALL ZGEMM( 'C', 'N', nbase, nbase, kdim, ONE, &
psi, kdmx, hpsi, kdmx, ZERO, hc, nvecx )
!
CALL reduce( 2 * nbase * nvecx, hc )
CALL reduce( 2*nbase*nvecx, hc )
!
IF ( lrot ) THEN
!
DO n = 1, nbase
!
e(n) = hc(n,n)
e(n) = REAL( hc(n,n) )
!
vc(n,n) = ONE
!
END DO
@ -276,7 +275,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
DO np = 1, notcnv
!
psi(:,:,nbase+np) = - ew(nbase+np) * psi(:,:,nbase+np)
psi(:,:,nbase+np) = - ew(nbase+np)*psi(:,:,nbase+np)
!
END DO
!
@ -335,9 +334,10 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
CALL h_psi_nc( ndmx, ndim, notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
!
IF ( lelfield ) &
CALL h_epsi_her( ndmx, ndim, notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
CALL h_epsi_her( ndmx, ndim, &
notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL s_psi_nc( ndmx, ndim, notcnv, psi(1,1,nb1), spsi(1,1,nb1) )
!
@ -355,9 +355,10 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
CALL h_psi( ndmx, ndim, notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
!
IF ( lelfield ) &
CALL h_epsi_her( ndmx, ndim, notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
CALL h_epsi_her( ndmx, ndim, &
notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL s_psi( ndmx, ndim, notcnv, psi(1,1,nb1), spsi(1,1,nb1) )
!
@ -368,20 +369,20 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
END IF
!
CALL cgramg1( kdmx, nvecx, kdim, &
nb1, nbase+notcnv, psi, spsi, hpsi )
nb1, nbase+notcnv, psi, spsi, hpsi )
!
END IF
!
! ... update the reduced hamiltonian
!
CALL start_clock( 'overlap' )
CALL start_clock( 'uspp' )
!
CALL ZGEMM( 'C', 'N', nbase+notcnv, notcnv, kdim, ONE, psi, &
kdmx, hpsi(1,1,nb1), kdmx, ZERO, hc(1,nb1), nvecx )
!
CALL reduce( 2 * nvecx * notcnv, hc(1,nb1) )
CALL reduce( 2*nvecx*notcnv, hc(1,nb1) )
!
CALL stop_clock( 'overlap' )
CALL stop_clock( 'uspp' )
!
nbase = nbase + notcnv
!
@ -389,7 +390,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
! ... the diagonal of hc must be strictly real
!
hc(n,n) = CMPLX( DBLE( hc(n,n) ), 0.D0 )
hc(n,n) = CMPLX( REAL( hc(n,n) ), 0.D0 )
!
DO m = n + 1, nbase
!
@ -446,7 +447,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
! ... last iteration, some roots not converged: return
!
WRITE( stdout, '(5X,"WARNING: ",I5, &
& " eigenvalues not converged")' ) notcnv
& " eigenvalues not converged")' ) notcnv
!
CALL stop_clock( 'last' )
!
@ -477,7 +478,8 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
DO n = 1, nbase
!
hc(n,n) = e(n)
hc(n,n) = REAL( e(n) )
!
vc(n,n) = ONE
!
END DO
@ -504,9 +506,11 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
IMPLICIT NONE
!
ALLOCATE( psi( ndmx, npol, nvecx ) )
ALLOCATE( hpsi( ndmx, npol, nvecx ) )
IF ( overlap ) ALLOCATE( spsi( ndmx, npol, nvecx ) )
ALLOCATE( psi( ndmx, npol, nvecx ) )
ALLOCATE( hpsi( ndmx, npol, nvecx ) )
!
IF ( uspp ) ALLOCATE( spsi( ndmx, npol, nvecx ) )
!
ALLOCATE( sc( nvecx, nvecx ) )
ALLOCATE( hc( nvecx, nvecx ) )
ALLOCATE( vc( nvecx, nvecx ) )
@ -517,9 +521,10 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
nbase = nvec
conv = .FALSE.
!
IF ( overlap ) spsi = ZERO
psi = ZERO
IF ( uspp ) spsi = ZERO
!
hpsi = ZERO
psi = ZERO
psi(:,:,1:nvec) = evc(:,:,1:nvec)
!
! ... hpsi contains h times the basis vectors
@ -529,14 +534,14 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
CALL h_psi_nc( ndmx, ndim, nvec, psi, hpsi )
!
IF ( lelfield ) CALL h_epsi_her( ndmx, ndim, nvec, psi, hpsi )
IF ( overlap ) CALL s_psi_nc( ndmx, ndim, nvec, psi, spsi )
IF ( uspp ) CALL s_psi_nc( ndmx, ndim, nvec, psi, spsi )
!
ELSE
!
CALL h_psi( ndmx, ndim, nvec, psi, hpsi )
!
IF ( lelfield ) CALL h_epsi_her( ndmx, ndim, nvec, psi, hpsi )
IF ( overlap ) CALL s_psi( ndmx, ndim, nvec, psi, spsi )
IF ( uspp ) CALL s_psi( ndmx, ndim, nvec, psi, spsi )
!
END IF
!
@ -552,7 +557,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
CALL reduce( 2*nbase*nvecx, hc )
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL ZGEMM( 'C', 'N', nbase, nbase, kdim, ONE, &
psi, kdmx, spsi, kdmx, ZERO, sc, nvecx )
@ -570,7 +575,8 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
DO n = 1, nbase
!
e(n) = hc(n,n)
e(n) = REAL( hc(n,n) )
!
vc(n,n) = ONE
!
END DO
@ -621,7 +627,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
! ... expand the basis set with new basis vectors ( H - e*S )|psi> ...
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL ZGEMM( 'N', 'N', kdim, notcnv, nbase, ONE, spsi, &
kdmx, vc, nvecx, ZERO, psi(1,1,nb1), kdmx )
@ -635,7 +641,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
DO np = 1, notcnv
!
psi(:,:,nbase+np) = - ew(nbase+np) * psi(:,:,nbase+np)
psi(:,:,nbase+np) = - ew(nbase+np)*psi(:,:,nbase+np)
!
END DO
!
@ -648,8 +654,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
IF ( noncolin ) THEN
!
CALL g_psi_nc( ndmx, ndim, notcnv, &
npol, psi(1,1,nb1), ew(nb1) )
CALL g_psi_nc( ndmx, ndim, notcnv, npol, psi(1,1,nb1), ew(nb1) )
!
ELSE
!
@ -665,7 +670,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
DO n = 1, notcnv
!
nbn = nbase+n
nbn = nbase + n
!
IF ( npol == 1 ) THEN
!
@ -695,36 +700,35 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
CALL h_psi_nc( ndmx, ndim, notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
!
IF ( lelfield ) &
CALL h_epsi_her( ndmx, ndim, &
notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
IF ( overlap ) &
CALL s_psi_nc( ndmx, ndim, notcnv, &
psi(1,1,nb1), spsi(1,1,nb1) )
CALL h_epsi_her( ndmx, ndim, &
notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
!
IF ( uspp ) &
CALL s_psi_nc( ndmx, ndim, notcnv, psi(1,1,nb1), spsi(1,1,nb1) )
!
ELSE
!
CALL h_psi( ndmx, ndim, notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
!
IF ( lelfield ) &
CALL h_epsi_her( ndmx, ndim, &
notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
IF ( overlap ) &
CALL s_psi( ndmx, ndim, notcnv, psi(1,1,nb1), spsi(1,1,nb1) )
CALL h_epsi_her( ndmx, ndim, &
notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
!
IF ( uspp ) &
CALL s_psi( ndmx, ndim, notcnv, psi(1,1,nb1), spsi(1,1,nb1) )
!
END IF
!
! ... update the reduced hamiltonian
!
CALL start_clock( 'overlap' )
CALL start_clock( 'uspp' )
!
CALL ZGEMM( 'C', 'N', nbase+notcnv, notcnv, kdim, ONE, psi, &
kdmx, hpsi(1,1,nb1), kdmx, ZERO, hc(1,nb1), nvecx )
!
CALL reduce( 2*nvecx*notcnv, hc(1,nb1) )
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL ZGEMM( 'C', 'N', nbase+notcnv, notcnv, kdim, ONE, psi, &
kdmx, spsi(1,1,nb1), kdmx, ZERO, sc(1,nb1), nvecx )
@ -738,7 +742,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
CALL reduce( 2*nvecx*notcnv, sc(1,nb1) )
!
CALL stop_clock( 'overlap' )
CALL stop_clock( 'uspp' )
!
nbase = nbase + notcnv
!
@ -746,8 +750,8 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
! ... the diagonal of hc and sc must be strictly real
!
hc(n,n) = CMPLX( DBLE( hc(n,n) ), 0.D0 )
sc(n,n) = CMPLX( DBLE( sc(n,n) ), 0.D0 )
hc(n,n) = CMPLX( REAL( hc(n,n) ), 0.D0 )
sc(n,n) = CMPLX( REAL( sc(n,n) ), 0.D0 )
!
DO m = n + 1, nbase
!
@ -817,7 +821,7 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
psi(:,:,1:nvec) = evc(:,:,1:nvec)
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL ZGEMM( 'N', 'N', kdim, nvec, nbase, ONE, spsi, &
kdmx, vc, nvecx, ZERO, psi(1,1,nvec+1), kdmx )
@ -841,7 +845,8 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
DO n = 1, nbase
!
hc(n,n) = e(n)
hc(n,n) = REAL( e(n) )
!
sc(n,n) = ONE
vc(n,n) = ONE
!
@ -858,7 +863,9 @@ SUBROUTINE cegterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
DEALLOCATE( vc )
DEALLOCATE( hc )
DEALLOCATE( sc )
IF ( overlap ) DEALLOCATE( spsi )
!
IF ( uspp ) DEALLOCATE( spsi )
!
DEALLOCATE( hpsi )
DEALLOCATE( psi )
!

68
PW/cgramg1.f90
View File

@ -26,28 +26,20 @@ SUBROUTINE cgramg1( lda, nvecx, n, start, finish, psi, spsi, hpsi )
!
IMPLICIT NONE
!
! ... first the dummy variables
!
INTEGER :: lda, n, nvecx, start, finish
INTEGER, INTENT(IN) :: &
lda, n, nvecx, start, finish
! input: leading dimension of the vectors
! input: physical dimension
! input: dimension of psi
! input: first vector to orthogonalize
! input: last vector to orthogonalize
COMPLEX(DP) :: psi(lda,nvecx), spsi(lda,nvecx), hpsi(lda,nvecx)
COMPLEX(DP), INTENT(INOUT) :: &
psi(lda,nvecx), spsi(lda,nvecx), hpsi(lda,nvecx)
! input/output: the vectors to be orthogonalized
!
! ... parameters
!
INTEGER, PARAMETER :: ierrx = 3
! maximum number of errors
!
! ... here the local variables
!
INTEGER :: vec, vecp, ierr
INTEGER :: vec, vecp
! counter on vectors
! counter on vectors
! counter on errors
REAL(DP) :: psi_norm
! the norm of a vector
REAL(DP), EXTERNAL :: DDOT
@ -85,17 +77,15 @@ SUBROUTINE cgramg1( lda, nvecx, n, start, finish, psi, spsi, hpsi )
!
ALLOCATE( ps( finish ) )
!
ierr = 0
!
DO vec = start, finish
!
IF ( vec > 1 ) THEN
!
CALL DGEMV( 'T', 2*n, vec-1, 2.D0, &
psi, 2*n, spsi(1,vec), 1, 0.D0, ps, 1 )
psi, 2*lda, spsi(1,vec), 1, 0.D0, ps, 1 )
!
IF ( gstart == 2 ) &
ps(1:vec-1) = ps(1:vec-1) - psi(1,1:vec-1) * spsi(1,vec)
ps(1:vec-1) = ps(1:vec-1) - psi(1,1:vec-1)*spsi(1,vec)
!
CALL reduce( ( vec - 1 ), ps )
!
@ -109,17 +99,17 @@ SUBROUTINE cgramg1( lda, nvecx, n, start, finish, psi, spsi, hpsi )
!
END IF
!
psi_norm = 2.D0 * DDOT( 2 * n, psi(1,vec), 1, spsi(1,vec), 1 )
psi_norm = 2.D0*DDOT( 2*n, psi(1,vec), 1, spsi(1,vec), 1 )
!
IF ( gstart == 2 ) psi_norm = psi_norm - psi(1,vec) * spsi(1,vec)
!
CALL reduce( 1, psi_norm )
!
IF ( psi_norm < 0.D0 ) THEN
IF ( psi_norm < eps8 ) THEN
!
WRITE( stdout, '(/,5X,"norm = ",F16.10,I4,/)' ) psi_norm, vec
!
CALL errore( 'cgramg1_gamma', 'negative norm in S ', 1 )
CALL errore( 'cgramg1_gamma', 'negative or zero norm in S ', 1 )
!
END IF
!
@ -129,16 +119,6 @@ SUBROUTINE cgramg1( lda, nvecx, n, start, finish, psi, spsi, hpsi )
hpsi(:,vec) = psi_norm * hpsi(:,vec)
spsi(:,vec) = psi_norm * spsi(:,vec)
!
IF ( psi_norm < eps8 ) THEN
!
ierr = ierr + 1
!
IF ( ierr <= ierrx ) CYCLE
!
CALL errore( 'cgramg1_gamma', 'absurd correction vector', vec )
!
END IF
!
END DO
!
DEALLOCATE( ps )
@ -155,27 +135,23 @@ SUBROUTINE cgramg1( lda, nvecx, n, start, finish, psi, spsi, hpsi )
!
COMPLEX(DP), ALLOCATABLE :: ps(:)
!
COMPLEX(DP) :: ZDOTC
!
!
ALLOCATE( ps( finish ) )
!
ierr = 0
!
DO vec = start, finish
!
IF ( vec > 1 ) THEN
!
CALL ZGEMV( 'C', n, vec-1, ONE, &
psi, n, spsi(1,vec), 1, ZERO, ps, 1 )
psi(1,1), lda, spsi(1,vec), 1, ZERO, ps, 1 )
!
CALL reduce( 2*( vec - 1 ), ps )
!
DO vecp = 1, ( vec - 1 )
DO vecp = 1, vec - 1
!
psi(:,vec) = psi(:,vec) - ps(vecp) * psi(:,vecp)
hpsi(:,vec) = hpsi(:,vec) - ps(vecp) * hpsi(:,vecp)
spsi(:,vec) = spsi(:,vec) - ps(vecp) * spsi(:,vecp)
psi(:,vec) = psi(:,vec) - ps(vecp)*psi(:,vecp)
hpsi(:,vec) = hpsi(:,vec) - ps(vecp)*hpsi(:,vecp)
spsi(:,vec) = spsi(:,vec) - ps(vecp)*spsi(:,vecp)
!
END DO
!
@ -185,11 +161,11 @@ SUBROUTINE cgramg1( lda, nvecx, n, start, finish, psi, spsi, hpsi )
!
CALL reduce( 1, psi_norm )
!
IF ( psi_norm < 0.D0 ) THEN
IF ( psi_norm < eps8 ) THEN
!
WRITE( stdout, '(/,5X,"norm = ",F16.10,I4,/)' ) psi_norm, vec
!
CALL errore( 'cgramg1_k', ' negative norm in S ', 1 )
CALL errore( 'cgramg1_k', 'negative or zero norm in S', 1 )
!
END IF
!
@ -199,16 +175,6 @@ SUBROUTINE cgramg1( lda, nvecx, n, start, finish, psi, spsi, hpsi )
hpsi(:,vec) = psi_norm * hpsi(:,vec)
spsi(:,vec) = psi_norm * spsi(:,vec)
!
IF ( psi_norm < eps8 ) THEN
!
ierr = ierr + 1
!
IF ( ierr <= ierrx ) CYCLE
!
CALL errore( 'cgramg1_k', ' absurd correction vector', vec )
!
END IF
!
END DO
!
DEALLOCATE( ps )

45
PW/regterg.f90
View File

@ -12,7 +12,7 @@
!
!----------------------------------------------------------------------------
SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
overlap, gstart, e, btype, notcnv, lrot, dav_iter )
uspp, gstart, e, btype, notcnv, lrot, dav_iter )
!----------------------------------------------------------------------------
!
! ... iterative solution of the eigenvalue problem:
@ -20,7 +20,7 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
! ... ( H - e S ) * evc = 0
!
! ... where H is an hermitean operator, e is a real scalar,
! ... S is an overlap matrix, evc is a complex vector
! ... S is an uspp matrix, evc is a complex vector
! ... (real wavefunctions with only half plane waves stored)
!
USE io_global, ONLY : stdout
@ -42,7 +42,7 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
REAL(DP), INTENT(IN) :: ethr
! energy threshold for convergence: root improvement is stopped,
! when two consecutive estimates of the root differ by less than ethr.
LOGICAL, INTENT(IN) :: overlap
LOGICAL, INTENT(IN) :: uspp
! if .FALSE. : S|psi> not needed
INTEGER, INTENT(IN) :: btype(nvec)
! band type ( 1 = occupied, 0 = empty )
@ -152,7 +152,7 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
CALL h_psi( ndmx, ndim, nvec, psi, hpsi )
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL s_psi( ndmx, ndim, nvec, psi, spsi )
!
@ -277,7 +277,7 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
CALL h_psi( ndmx, ndim, notcnv, psi(1,nb1), hpsi(1,nb1) )
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL s_psi( ndmx, ndim, notcnv, psi(1,nb1), spsi(1,nb1) )
!
@ -291,7 +291,7 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
! ... update the reduced hamiltonian
!
CALL start_clock( 'overlap' )
CALL start_clock( 'uspp' )
!
CALL DGEMM( 'T', 'N', nbase+notcnv, notcnv, ndim2, 2.D0, psi, &
ndmx2, hpsi(1,nb1), ndmx2, 0.D0, hr(1,nb1), nvecx )
@ -302,7 +302,7 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
CALL reduce( nvecx * notcnv, hr(1,nb1) )
!
CALL stop_clock( 'overlap' )
CALL stop_clock( 'uspp' )
!
nbase = nbase + notcnv
!
@ -421,9 +421,11 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
IMPLICIT NONE
!
ALLOCATE( psi( ndmx, nvecx ) )
ALLOCATE( psi( ndmx, nvecx ) )
ALLOCATE( hpsi( ndmx, nvecx ) )
IF ( overlap ) ALLOCATE( spsi( ndmx, nvecx ) )
!
IF ( uspp ) ALLOCATE( spsi( ndmx, nvecx ) )
!
ALLOCATE( sr( nvecx, nvecx ) )
ALLOCATE( hr( nvecx, nvecx ) )
ALLOCATE( vr( nvecx, nvecx ) )
@ -436,16 +438,17 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
nbase = nvec
conv = .FALSE.
!
IF ( overlap ) spsi = ZERO
psi = ZERO
IF ( uspp ) spsi = ZERO
!
hpsi = ZERO
psi = ZERO
psi(:,1:nvec) = evc(:,1:nvec)
!
! ... hpsi contains h times the basis vectors
!
CALL h_psi( ndmx, ndim, nvec, psi, hpsi )
!
IF ( overlap ) CALL s_psi( ndmx, ndim, nvec, psi, spsi )
IF ( uspp ) CALL s_psi( ndmx, ndim, nvec, psi, spsi )
!
! ... hr contains the projection of the hamiltonian onto the reduced
! ... space vr contains the eigenvectors of hr
@ -462,7 +465,7 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
CALL reduce( nbase*nvecx, hr )
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL DGEMM( 'T', 'N', nbase, nbase, ndim2, 2.D0, &
psi, ndmx2, spsi, ndmx2, 0.D0, sr, nvecx )
@ -537,7 +540,7 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
! ... expand the basis set with new basis vectors ( H - e*S )|psi> ...
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL DGEMM( 'N', 'N', ndim2, notcnv, nbase, 1.D0, &
spsi, ndmx2, vr, nvecx, 0.D0, psi(1,nb1), ndmx2 )
@ -590,12 +593,12 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
CALL h_psi( ndmx, ndim, notcnv, psi(1,nb1), hpsi(1,nb1) )
!
IF ( overlap ) &
IF ( uspp ) &
CALL s_psi( ndmx, ndim, notcnv, psi(1,nb1), spsi(1,nb1) )
!
! ... update the reduced hamiltonian
!
CALL start_clock( 'overlap' )
CALL start_clock( 'uspp' )
!
CALL DGEMM( 'T', 'N', nbase+notcnv, notcnv, ndim2, 2.D0, psi, &
ndmx2, hpsi(1,nb1), ndmx2, 0.D0, hr(1,nb1), nvecx )
@ -606,7 +609,7 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
CALL reduce( nvecx * notcnv, hr(1,nb1) )
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL DGEMM( 'T', 'N', nbase+notcnv, notcnv, ndim2, 2.D0, psi, &
ndmx2, spsi(1,nb1), ndmx2, 0.D0, sr(1,nb1), nvecx )
@ -628,7 +631,7 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
CALL reduce( nvecx*notcnv, sr(1,nb1) )
!
CALL stop_clock( 'overlap' )
CALL stop_clock( 'uspp' )
!
nbase = nbase + notcnv
!
@ -702,7 +705,7 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
!
psi(:,1:nvec) = evc(:,1:nvec)
!
IF ( overlap ) THEN
IF ( uspp ) THEN
!
CALL DGEMM( 'N', 'N', ndim2, nvec, nbase, 1.D0, spsi, &
ndmx2, vr, nvecx, 0.D0, psi(1,nvec+1), ndmx2 )
@ -743,7 +746,9 @@ SUBROUTINE regterg( ndim, ndmx, nvec, nvecx, evc, ethr, &
DEALLOCATE( vr )
DEALLOCATE( hr )
DEALLOCATE( sr )
IF ( overlap ) DEALLOCATE( spsi )
!
IF ( uspp ) DEALLOCATE( spsi )
!
DEALLOCATE( hpsi )
DEALLOCATE( psi )
!

23
PW/setup.f90
View File

@ -58,7 +58,7 @@ SUBROUTINE setup()
USE ktetra, ONLY : nk1, nk2, nk3, k1, k2, k3, &
tetra, ntetra, ltetra
USE symme, ONLY : s, t_rev, irt, ftau, nsym, invsym
USE atom, ONLY : r, oc, chi, nchi, lchi, jchi, mesh
USE atom, ONLY : oc, chi, nchi, lchi, jchi, mesh
USE pseud, ONLY : zp
USE wvfct, ONLY : nbnd, nbndx, gamma_only
USE control_flags, ONLY : tr2, ethr, alpha0, beta0, lscf, lmd, lpath, &
@ -84,22 +84,15 @@ SUBROUTINE setup()
!
IMPLICIT NONE
!
! ... local variables
!
REAL(DP), PARAMETER :: &
eps = 1.0D-12 ! small number
INTEGER :: &
na, &!
ir, &!
nt, &!
input_nks, &!
nrot, &!
iter, &!
ierr, &!
irot, &!
isym, &!
ipol, &!
jpol, &!
tipo, &!
is, &!
nb, &!
@ -111,13 +104,13 @@ SUBROUTINE setup()
minus_q, &!
ltest !
REAL(DP) :: &
vionl, & !
vionl, &!
iocc !
INTEGER, EXTERNAL :: &
n_atom_wfc, &!
set_Hubbard_l
LOGICAL, EXTERNAL :: &
lchk_tauxk ! tests that atomic coordinates do not overlap
lchk_tauxk ! tests that atomic coordinates do not overlap
!
!
ALLOCATE( m_loc( 3, nat ) )
@ -854,14 +847,14 @@ FUNCTION n_atom_wfc( nat, npsx, ityp, nchix, nchi, oc, lchi, jchi )
!
USE kinds, ONLY : DP
use noncollin_module, ONLY : noncolin
use spin_orb, ONLY : lspinorb, so
use spin_orb, ONLY : so
!
IMPLICIT NONE
!
INTEGER :: n_atom_wfc
INTEGER :: nat, npsx, ityp(nat), nchix, nchi(npsx), lchi(nchix,npsx)
INTEGER :: n_atom_wfc
INTEGER :: nat, npsx, ityp(nat), nchix, nchi(npsx), lchi(nchix,npsx)
REAL(DP) :: oc(nchix,npsx), jchi(nchix,npsx)
INTEGER :: na, nt, n
INTEGER :: na, nt, n
!
!
n_atom_wfc = 0
@ -1068,7 +1061,7 @@ SUBROUTINE check_para_diag_efficiency()
IF ( time_para < time_serial ) THEN
!
use_para_diago = .TRUE.
para_diago_dim = dim
!
EXIT
!