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

Phonon cleanup. The preconditioning variable eprec is calculated in a single 

place. h_diag is calculated outside the loop on the perturbations.


git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@5412 c92efa57-630b-4861-b058-cf58834340f0
This commit is contained in:
dalcorso 2009-02-14 07:46:55 +00:00
parent 94f4212854
commit 097596216d
8 changed files with 75 additions and 126 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

3
PH/allocate_phq.f90
View File

@ -36,7 +36,7 @@ subroutine allocate_phq
dpqq_so, alphasum, alphasum_nc, becsum_nc, &
becp1, becp1_nc, alphap, alphap_nc
USE efield_mod, ONLY : zstareu, zstareu0, zstarue0, zstarue0_rec, zstarue
USE eqv, ONLY : dpsi, evq, vlocq, dmuxc, dvpsi
USE eqv, ONLY : dpsi, evq, vlocq, dmuxc, dvpsi, eprec
USE units_ph, ONLY : this_pcxpsi_is_on_file, this_dvkb3_is_on_file
USE dynmat, ONLY : dyn00, dyn, dyn_rec, w2
USE modes, ONLY : u, ubar, rtau, max_irr_dim, npert, t, tmq
@ -67,6 +67,7 @@ subroutine allocate_phq
!
allocate (vlocq ( ngm , ntyp))
allocate (dmuxc ( nrxx , nspin , nspin))
allocate (eprec ( nbnd, nksq) )
!
allocate (eigqts ( nat))
allocate (rtau ( 3, 48, nat))

3
PH/deallocate_phq.f90
View File

@ -28,7 +28,7 @@ subroutine deallocate_phq
vsgga, segni
USE gamma_gamma, ONLY : with_symmetry, has_equivalent, equiv_atoms, &
n_equiv_atoms
USE eqv, ONLY : dmuxc, vlocq, dpsi, dvpsi, evq
USE eqv, ONLY : dmuxc, vlocq, dpsi, dvpsi, evq, eprec
USE nlcc_ph, ONLY : drc
USE units_ph, ONLY : this_dvkb3_is_on_file, this_pcxpsi_is_on_file
USE dynmat, ONLY : dyn00, dyn_rec, dyn, w2
@ -51,6 +51,7 @@ subroutine deallocate_phq
!
if(allocated(vlocq)) deallocate (vlocq)
if(allocated(dmuxc)) deallocate (dmuxc)
if(allocated(eprec)) deallocate (eprec)
!
if(allocated(ikks)) deallocate (ikks)
if(allocated(ikqs)) deallocate (ikqs)

26
PH/dvpsi_e.f90
View File

@ -32,7 +32,7 @@ subroutine dvpsi_e (ik, ipol)
USE uspp, ONLY : okvan, nkb, vkb, qq, qq_so, deeq, deeq_nc
USE uspp_param, ONLY : nh
USE ramanm, ONLY : eth_rps
USE eqv, ONLY : dpsi, dvpsi
USE eqv, ONLY : dpsi, dvpsi, eprec
USE phus, ONLY : becp1, becp1_nc
USE qpoint, ONLY : npwq, nksq
USE units_ph, ONLY : this_pcxpsi_is_on_file, lrcom, iucom, &
@ -51,7 +51,7 @@ subroutine dvpsi_e (ik, ipol)
integer :: ig, na, ibnd, jbnd, ikb, jkb, nt, lter, ih, jh, ijkb0, nrec
! counters
real(DP), allocatable :: gk (:,:), h_diag (:,:), eprec (:)
real(DP), allocatable :: gk (:,:), h_diag (:,:)
! the derivative of |k+G|
real(DP) :: anorm, thresh
! preconditioning cut-off
@ -82,7 +82,6 @@ subroutine dvpsi_e (ik, ipol)
allocate (aux ( npwx*npol ))
allocate (gk ( 3, npwx))
allocate (h_diag( npwx*npol, nbnd))
allocate (eprec( nbnd))
if (nkb > 0) then
IF (noncolin) THEN
allocate (becp2_nc (nkb, npol, nbnd))
@ -305,30 +304,14 @@ subroutine dvpsi_e (ik, ipol)
! Now solve the linear systems (H-e_vS)*P_c(x*psi_v)=P_c^+ [H-e_vS,x]*psi_v
!
thresh = eth_rps
do ibnd = 1, nbnd_occ (ik)
conv_root = .true.
aux=(0.d0,0.d0)
do ig = 1, npwq
aux(ig) = g2kin (ig) * evc (ig, ibnd)
enddo
IF (noncolin) THEN
do ig = 1, npwq
aux(ig+npwx) = g2kin (ig) * evc (ig+npwx, ibnd)
enddo
END IF
eprec (ibnd) = 1.35d0 * ZDOTC (npwx*npol, evc (1, ibnd), 1, aux, 1)
enddo
#ifdef __PARA
call mp_sum ( eprec( 1:nbnd_occ(ik) ), intra_pool_comm )
#endif
h_diag=0.d0
do ibnd = 1, nbnd_occ (ik)
do ig = 1, npwq
h_diag (ig, ibnd) = 1.d0 / max (1.0d0, g2kin (ig) / eprec (ibnd) )
h_diag (ig, ibnd) = 1.d0 / max (1.0d0, g2kin (ig) / eprec (ibnd,ik) )
enddo
IF (noncolin) THEN
do ig = 1, npwq
h_diag (ig+npwx, ibnd) = 1.d0/max(1.0d0,g2kin(ig)/eprec(ibnd))
h_diag (ig+npwx, ibnd) = 1.d0/max(1.0d0,g2kin(ig)/eprec(ibnd,ik))
enddo
END IF
enddo
@ -388,7 +371,6 @@ subroutine dvpsi_e (ik, ipol)
ENDIF
END IF
deallocate (eprec)
deallocate (h_diag)
deallocate (work)
deallocate (aux)

18
PH/pcgreen.f90
View File

@ -21,7 +21,7 @@ subroutine pcgreen (avg_iter, thresh, ik, et_ )
USE wavefunctions_module, ONLY: evc
USE mp_global, ONLY: intra_pool_comm
USE mp, ONLY: mp_sum
USE eqv, ONLY: dpsi, dvpsi
USE eqv, ONLY: dpsi, dvpsi, eprec
USE control_ph, ONLY : nbnd_occ
implicit none
@ -50,9 +50,8 @@ subroutine pcgreen (avg_iter, thresh, ik, et_ )
real(DP) :: anorm
! the norm of the error
real(DP) , allocatable :: h_diag(:,:), eprec(:)
real(DP) , allocatable :: h_diag(:,:)
! the diagonal part of the Hamiltonian
! cut-off for preconditioning
complex(DP) :: ZDOTC
! the scalar product function
@ -64,7 +63,6 @@ subroutine pcgreen (avg_iter, thresh, ik, et_ )
allocate (h_diag ( npwx, nbnd ))
allocate (auxg ( npwx ))
allocate (eprec ( nbnd ))
allocate (ps ( nbnd, nbnd ))
!
@ -98,16 +96,7 @@ subroutine pcgreen (avg_iter, thresh, ik, et_ )
!
do ibnd = 1, nbnd_occ (ik)
do ig = 1, npw
auxg (ig) = g2kin (ig) * evc(ig, ibnd)
enddo
eprec (ibnd) = 1.35d0 * ZDOTC (npw, evc(1, ibnd), 1, auxg, 1)
enddo
#ifdef __PARA
call mp_sum ( eprec( 1:nbnd_occ(ik) ), intra_pool_comm )
#endif
do ibnd = 1, nbnd_occ (ik)
do ig = 1, npw
h_diag (ig, ibnd) = 1.d0 / max (1.0d0, g2kin (ig) / eprec (ibnd) )
h_diag (ig, ibnd) = 1.d0 / max (1.0d0, g2kin (ig) / eprec (ibnd,ik) )
enddo
enddo
conv_root = .true.
@ -122,7 +111,6 @@ subroutine pcgreen (avg_iter, thresh, ik, et_ )
ik,ibnd,anorm
deallocate (ps)
deallocate (eprec)
deallocate (auxg)
deallocate (h_diag)

42
PH/phq_init.f90
View File

@ -51,21 +51,23 @@ SUBROUTINE phq_init()
USE noncollin_module, ONLY : noncolin, npol
USE uspp, ONLY : okvan, vkb
USE uspp_param, ONLY : upf
USE eqv, ONLY : vlocq, evq
USE eqv, ONLY : vlocq, evq, eprec
USE phus, ONLY : becp1, becp1_nc, alphap, alphap_nc, dpqq, &
dpqq_so
USE nlcc_ph, ONLY : nlcc_any
USE control_ph, ONLY : zue, epsil, lgamma, all_done
USE control_ph, ONLY : zue, epsil, lgamma, all_done, nbnd_occ
USE units_ph, ONLY : lrwfc, iuwfc
USE qpoint, ONLY : xq, igkq, npwq, nksq, eigqts, ikks, ikqs
USE paw_onecenter, ONLY : PAW_potential, PAW_symmetrize
USE paw_variables, ONLY : okpaw, ddd_paw
USE mp_global, ONLY : intra_pool_comm
USE mp, ONLY : mp_sum
!
IMPLICIT NONE
!
! ... local variables
!
INTEGER :: nt, ik, ikq, ipol, ibnd, ikk, na, ig
INTEGER :: nt, ik, ikq, ipol, ibnd, ikk, na, ig, irr, imode0
! counter on atom types
! counter on k points
! counter on k+q points
@ -78,6 +80,7 @@ SUBROUTINE phq_init()
! the argument of the phase
COMPLEX(DP), ALLOCATABLE :: aux1(:,:)
! used to compute alphap
COMPLEX(DP), EXTERNAL :: ZDOTC
!
!
IF (all_done) RETURN
@ -203,13 +206,36 @@ SUBROUTINE phq_init()
END IF
END DO
!
! ... if there is only one k-point the k+q wavefunctions are
! ... read once here
!
IF ( nksq == 1 .AND. .NOT. lgamma ) &
IF ( .NOT. lgamma ) &
CALL davcio( evq, lrwfc, iuwfc, ikq, -1 )
!
! diagonal elements of the unperturbed Hamiltonian,
! needed for preconditioning
!
do ig = 1, npwq
g2kin (ig) = ( (xk (1,ikq) + g (1, igkq(ig)) ) **2 + &
(xk (2,ikq) + g (2, igkq(ig)) ) **2 + &
(xk (3,ikq) + g (3, igkq(ig)) ) **2 ) * tpiba2
enddo
aux1=(0.d0,0.d0)
DO ig = 1, npwq
aux1 (ig,1:nbnd_occ(ikk)) = g2kin (ig) * evq (ig, 1:nbnd_occ(ikk))
END DO
IF (noncolin) THEN
DO ig = 1, npwq
aux1 (ig+npwx,1:nbnd_occ(ikk)) = g2kin (ig)* &
evq (ig+npwx, 1:nbnd_occ(ikk))
END DO
END IF
DO ibnd=1,nbnd_occ(ikk)
eprec (ibnd,ik) = 1.35d0 * ZDOTC(npwx*npol,evq(1,ibnd),1,aux1(1,ibnd),1)
END DO
!
END DO
#ifdef __PARA
CALL mp_sum ( eprec, intra_pool_comm )
#endif
!
DEALLOCATE( aux1 )
!

44
PH/solve_e.f90
View File

@ -42,7 +42,7 @@ subroutine solve_e
USE paw_variables, ONLY : okpaw
USE paw_onecenter, ONLY : paw_dpotential, paw_desymmetrize
USE eqv, ONLY : dpsi, dvpsi
USE eqv, ONLY : dpsi, dvpsi, eprec
USE modes, ONLY : max_irr_dim
USE units_ph, ONLY : lrdwf, iudwf, lrwfc, iuwfc, lrdrho, &
iudrho, iunrec, this_pcxpsi_is_on_file
@ -63,9 +63,8 @@ subroutine solve_e
! anorm : the norm of the error
! averlt: average number of iterations
! dr2 : self-consistency error
real(DP), allocatable :: h_diag (:,:), eprec(:)
real(DP), allocatable :: h_diag (:,:)
! h_diag: diagonal part of the Hamiltonian
! eprec : array fo preconditioning
complex(DP) , allocatable, target :: &
dvscfin (:,:,:) ! change of the scf potential (input)
@ -110,7 +109,6 @@ subroutine solve_e
allocate (ps (nbnd,nbnd))
ps (:,:) = (0.d0, 0.d0)
allocate (h_diag(npwx*npol, nbnd))
allocate (eprec(nbnd))
if (rec_code == -20.and.recover) then
! restarting in Electric field calculation
CALL read_rec(dr2, iter0, dvscfin, dvscfins, 3)
@ -172,6 +170,17 @@ subroutine solve_e
(xk (2,ik ) + g (2,igkq (ig)) ) **2 + &
(xk (3,ik ) + g (3,igkq (ig)) ) **2 ) * tpiba2
enddo
h_diag=0.d0
do ibnd = 1, nbnd_occ (ik)
do ig = 1, npw
h_diag(ig,ibnd)=1.d0/max(1.0d0,g2kin(ig)/eprec(ibnd,ik))
enddo
IF (noncolin) THEN
do ig = 1, npw
h_diag(ig+npwx,ibnd)=1.d0/max(1.0d0,g2kin(ig)/eprec(ibnd,ik))
enddo
END IF
enddo
!
do ipol = 1, 3
!
@ -257,32 +266,6 @@ subroutine solve_e
! iterative solution of the linear system (H-e)*dpsi=dvpsi
! dvpsi=-P_c+ (dvbare+dvscf)*psi , dvscf fixed.
!
do ibnd = 1, nbnd_occ (ik)
auxg=(0.d0,0.d0)
do ig = 1, npw
auxg (ig) = g2kin (ig) * evc (ig, ibnd)
enddo
IF (noncolin) THEN
do ig = 1, npw
auxg (ig+npwx) = g2kin (ig) * evc (ig+npwx, ibnd)
enddo
END if
eprec (ibnd) = 1.35d0*ZDOTC(npwx*npol,evc(1,ibnd),1,auxg,1)
enddo
#ifdef __PARA
call mp_sum ( eprec( 1:nbnd_occ(ik) ), intra_pool_comm )
#endif
h_diag=0.d0
do ibnd = 1, nbnd_occ (ik)
do ig = 1, npw
h_diag(ig,ibnd)=1.d0/max(1.0d0,g2kin(ig)/eprec(ibnd))
enddo
IF (noncolin) THEN
do ig = 1, npw
h_diag(ig+npwx,ibnd)=1.d0/max(1.0d0,g2kin(ig)/eprec(ibnd))
enddo
END IF
enddo
conv_root = .true.
@ -415,7 +398,6 @@ subroutine solve_e
enddo
155 continue
deallocate (eprec)
deallocate (h_diag)
deallocate (ps)
deallocate (aux1)

20
PH/solve_e_fpol.f90
View File

@ -37,7 +37,7 @@ subroutine solve_e_fpol ( iw )
USE wvfct, ONLY : npw, npwx, nbnd, igk, g2kin, et
USE uspp, ONLY : okvan, vkb
USE uspp_param, ONLY : nhm
USE eqv, ONLY : dpsi, dvpsi
USE eqv, ONLY : dpsi, dvpsi, eprec
USE control_ph, ONLY : nmix_ph, tr2_ph, alpha_mix, convt, &
nbnd_occ, lgamma, niter_ph, &
rec_code, flmixdpot
@ -60,9 +60,6 @@ subroutine solve_e_fpol ( iw )
! the eigenvalues plus imaginary frequency
! the diagonal part of the Hamiltonian which becomes complex now
real(DP), allocatable :: eprec(:)
! real(DP), allocatable :: h_diag (:,:), eprec(:)
! eprec : array fo preconditioning
complex(DP) , allocatable, target :: &
dvscfin (:,:,:) ! change of the scf potential (input)
@ -106,7 +103,6 @@ subroutine solve_e_fpol ( iw )
allocate (ps (nbnd,nbnd))
ps (:,:) = (0.d0, 0.d0)
allocate (h_diag(npwx, nbnd))
allocate (eprec(nbnd))
allocate (etc(nbnd, nkstot))
etc(:,:) = CMPLX( et(:,:), iw )
@ -254,23 +250,14 @@ subroutine solve_e_fpol ( iw )
! iterative solution of the linear system (H-e)*dpsi=dvpsi
! dvpsi=-P_c+ (dvbare+dvscf)*psi , dvscf fixed.
!
do ibnd = 1, nbnd_occ (ik)
do ig = 1, npw
auxg (ig) = g2kin (ig) * evc (ig, ibnd)
enddo
eprec (ibnd) = 1.35d0*ZDOTC(npwq,evc(1,ibnd),1,auxg,1)
enddo
#ifdef __PARA
call mp_sum ( eprec( 1:nbnd_occ(ik) ), intra_pool_comm )
#endif
do ibnd = 1, nbnd_occ (ik)
!
if ( (abs(iw).lt.0.05) .or. (abs(iw).gt.1.d0) ) then
!
do ig = 1, npw
! h_diag(ig,ibnd)=1.d0/max(1.0d0,g2kin(ig)/eprec(ibnd))
! h_diag(ig,ibnd)=1.d0/max(1.0d0,g2kin(ig)/eprec(ibnd,ik))
h_diag(ig,ibnd)=CMPLX(1.d0, 0.d0) / &
CMPLX( max(1.0d0,g2kin(ig)/eprec(ibnd))-et(ibnd,ik),-iw )
CMPLX( max(1.0d0,g2kin(ig)/eprec(ibnd,ik))-et(ibnd,ik),-iw )
end do
else
do ig = 1, npw
@ -402,7 +389,6 @@ subroutine solve_e_fpol ( iw )
if (convt) goto 155
enddo
155 continue
deallocate (eprec)
deallocate (h_diag)
deallocate (ps)
deallocate (aux1)

45
PH/solve_linter.f90
View File

@ -55,7 +55,7 @@ subroutine solve_linter (irr, imode0, npe, drhoscf)
this_pcxpsi_is_on_file
USE output, ONLY : fildrho, fildvscf
USE phus, ONLY : int3_paw, becsumort
USE eqv, ONLY : dvpsi, dpsi, evq
USE eqv, ONLY : dvpsi, dpsi, evq, eprec
USE qpoint, ONLY : xq, npwq, igkq, nksq, ikks, ikqs
USE modes, ONLY : npert, u, t, max_irr_dim, irotmq, tmq, &
minus_q, irgq, nsymq, rtau
@ -73,9 +73,8 @@ subroutine solve_linter (irr, imode0, npe, drhoscf)
complex(DP) :: drhoscf (nrxx, nspin, npe)
! output: the change of the scf charge
real(DP) , allocatable :: h_diag (:,:),eprec (:)
real(DP) , allocatable :: h_diag (:,:)
! h_diag: diagonal part of the Hamiltonian
! eprec : array for preconditioning
real(DP) :: thresh, anorm, averlt, dr2
! thresh: convergence threshold
! anorm : the norm of the error
@ -153,7 +152,6 @@ subroutine solve_linter (irr, imode0, npe, drhoscf)
IF (noncolin) allocate (dbecsum_nc (nhm,nhm, nat , nspin , npe))
allocate (aux1 ( nrxxs, npol))
allocate (h_diag ( npwx*npol, nbnd))
allocate (eprec ( nbnd))
!
if (rec_code > 2.and.recover) then
! restart from Phonon calculation
@ -240,6 +238,18 @@ subroutine solve_linter (irr, imode0, npe, drhoscf)
(xk (2,ikq) + g (2, igkq(ig)) ) **2 + &
(xk (3,ikq) + g (3, igkq(ig)) ) **2 ) * tpiba2
enddo
h_diag=0.d0
do ibnd = 1, nbnd_occ (ikk)
do ig = 1, npwq
h_diag(ig,ibnd)=1.d0/max(1.0d0,g2kin(ig)/eprec(ibnd,ik))
enddo
IF (noncolin) THEN
do ig = 1, npwq
h_diag(ig+npwx,ibnd)=1.d0/max(1.0d0,g2kin(ig)/eprec(ibnd,ik))
enddo
END IF
enddo
!
! diagonal elements of the unperturbed hamiltonian
!
@ -396,32 +406,6 @@ subroutine solve_linter (irr, imode0, npe, drhoscf)
! iterative solution of the linear system (H-eS)*dpsi=dvpsi,
! dvpsi=-P_c^+ (dvbare+dvscf)*psi , dvscf fixed.
!
do ibnd = 1, nbnd_occ (ikk)
auxg=(0.d0,0.d0)
do ig = 1, npwq
auxg (ig) = g2kin (ig) * evq (ig, ibnd)
enddo
IF (noncolin) THEN
do ig = 1, npwq
auxg (ig+npwx) = g2kin (ig) * evq (ig+npwx, ibnd)
enddo
END IF
eprec (ibnd) = 1.35d0 * ZDOTC (npwx*npol,evq(1,ibnd),1,auxg, 1)
enddo
#ifdef __PARA
call mp_sum ( eprec( 1:nbnd_occ (ikk) ), intra_pool_comm )
#endif
h_diag=0.d0
do ibnd = 1, nbnd_occ (ikk)
do ig = 1, npwq
h_diag(ig,ibnd)=1.d0/max(1.0d0,g2kin(ig)/eprec(ibnd))
enddo
IF (noncolin) THEN
do ig = 1, npwq
h_diag(ig+npwx,ibnd)=1.d0/max(1.0d0,g2kin(ig)/eprec(ibnd))
enddo
END IF
enddo
conv_root = .true.
call cgsolve_all (ch_psi_all, cg_psi, et(1,ikk), dvpsi, dpsi, &
@ -646,7 +630,6 @@ subroutine solve_linter (irr, imode0, npe, drhoscf)
if (convt.and.nlcc_any) call addnlcc (imode0, drhoscfh, npe)
if (lmetq0) deallocate (ldoss)
if (lmetq0) deallocate (ldos)
deallocate (eprec)
deallocate (h_diag)
deallocate (aux1)
deallocate (dbecsum)