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

more scf_mod changes 

git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4398 c92efa57-630b-4861-b058-cf58834340f0
This commit is contained in:
degironc 2007-11-06 10:26:07 +00:00
parent 8c462e6c22
commit 79f0fd08c0
12 changed files with 164 additions and 127 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
Modules/make.depend
View File

@ -89,6 +89,7 @@ ions_nose.o : ions_base.o
ions_nose.o : kind.o
ions_nose.o : parameters.o
ions_nose.o : timestep.o
kplusg.o : kind.o
metadyn_base.o : basic_algebra_routines.o
metadyn_base.o : cell_base.o
metadyn_base.o : constants.o
@ -255,6 +256,7 @@ uspp.o : kind.o
uspp.o : parameters.o
uspp.o : pseudo_types.o
uspp.o : random_numbers.o
vms.o : kind.o
vxc_t.o : functionals.o
vxc_t.o : io_global.o
vxc_t.o : kind.o

29
PH/make.depend
View File

@ -109,6 +109,7 @@ bcast_ph_input.o : ramanm.o
bcast_ph_input1.o : ../Modules/mp.o
bcast_ph_input1.o : ../PW/pwcom.o
bcast_ph_input1.o : phcom.o
bicgstab.o : ../Modules/kind.o
ccg_psi.o : ../Modules/kind.o
cch_psi_all.o : ../Modules/kind.o
cch_psi_all.o : ../PW/becmod.o
@ -350,6 +351,12 @@ dynmatrix.o : ../PW/noncol.o
dynmatrix.o : ../PW/pwcom.o
dynmatrix.o : phcom.o
dynmatrix.o : ramanm.o
ec_rpa.o : ../Modules/io_files.o
ec_rpa.o : ../Modules/io_global.o
ec_rpa.o : ../Modules/kind.o
ec_rpa.o : ../Modules/random_numbers.o
ec_rpa.o : ../PW/pwcom.o
ec_rpa.o : phcom.o
ef_shift.o : ../Modules/cell_base.o
ef_shift.o : ../Modules/io_global.o
ef_shift.o : ../Modules/kind.o
@ -391,6 +398,7 @@ find_mode_sym.o : ../PW/noncol.o
find_mode_sym.o : ../PW/pwcom.o
find_mode_sym.o : phcom.o
generate_dynamical_matrix_c.o : ../Modules/kind.o
gl_weight.o : phcom.o
gmressolve_all.o : ../Modules/kind.o
h_psiq.o : ../Modules/kind.o
h_psiq.o : ../Modules/wavefunctions.o
@ -480,6 +488,7 @@ phonon.o : phcom.o
phonon.o : ramanm.o
phq_init.o : ../Modules/atom.o
phq_init.o : ../Modules/constants.o
phq_init.o : ../Modules/grid_paw_variables.o
phq_init.o : ../Modules/io_files.o
phq_init.o : ../Modules/io_global.o
phq_init.o : ../Modules/ions_base.o
@ -488,6 +497,7 @@ phq_init.o : ../Modules/uspp.o
phq_init.o : ../Modules/wavefunctions.o
phq_init.o : ../PW/noncol.o
phq_init.o : ../PW/pwcom.o
phq_init.o : ../PW/rad_paw_routines.o
phq_init.o : phcom.o
phq_readin.o : ../Modules/constants.o
phq_readin.o : ../Modules/control_flags.o
@ -689,6 +699,17 @@ solve_e_nscf.o : ../Modules/wavefunctions.o
solve_e_nscf.o : ../PW/becmod.o
solve_e_nscf.o : ../PW/pwcom.o
solve_e_nscf.o : phcom.o
solve_edv.o : ../Modules/check_stop.o
solve_edv.o : ../Modules/control_flags.o
solve_edv.o : ../Modules/io_files.o
solve_edv.o : ../Modules/io_global.o
solve_edv.o : ../Modules/ions_base.o
solve_edv.o : ../Modules/kind.o
solve_edv.o : ../Modules/uspp.o
solve_edv.o : ../Modules/wavefunctions.o
solve_edv.o : ../PW/becmod.o
solve_edv.o : ../PW/pwcom.o
solve_edv.o : phcom.o
solve_linter.o : ../Modules/check_stop.o
solve_linter.o : ../Modules/constants.o
solve_linter.o : ../Modules/io_files.o
@ -782,6 +803,10 @@ transform_int_so.o : ../PW/pwcom.o
transform_int_so.o : phcom.o
trntnsc.o : ../Modules/kind.o
trntnsr_3.o : ../Modules/kind.o
vc_dv.o : ../Modules/functionals.o
vc_dv.o : ../Modules/kind.o
vc_dv.o : ../PW/pwcom.o
vc_dv.o : phcom.o
write_dyn_on_file.o : ../Modules/kind.o
write_epsilon_and_zeu.o : ../Modules/io_global.o
write_epsilon_and_zeu.o : ../Modules/kind.o
@ -829,6 +854,7 @@ allocate_part.o : ../include/f_defs.h
allocate_phq.o : ../include/f_defs.h
bcast_ph_input.o : ../include/f_defs.h
bcast_ph_input1.o : ../include/f_defs.h
bicgstab.o : ../include/f_defs.h
cch_psi_all.o : ../include/f_defs.h
cgsolve_all.o : ../include/f_defs.h
ch_psi_all.o : ../include/f_defs.h
@ -865,6 +891,7 @@ dynmat0.o : ../include/f_defs.h
dynmat_us.o : ../include/f_defs.h
dynmatcc.o : ../include/f_defs.h
dynmatrix.o : ../include/f_defs.h
ec_rpa.o : ../include/f_defs.h
ef_shift.o : ../include/f_defs.h
el_opt.o : ../include/f_defs.h
elphon.o : ../include/f_defs.h
@ -915,6 +942,7 @@ solve_e.o : ../include/f_defs.h
solve_e2.o : ../include/f_defs.h
solve_e_fpol.o : ../include/f_defs.h
solve_e_nscf.o : ../include/f_defs.h
solve_edv.o : ../include/f_defs.h
solve_linter.o : ../include/f_defs.h
star_q.o : ../include/f_defs.h
sym_and_write_zue.o : ../include/f_defs.h
@ -937,6 +965,7 @@ transform_int_nc.o : ../include/f_defs.h
transform_int_so.o : ../include/f_defs.h
trntnsc.o : ../include/f_defs.h
trntnsr_3.o : ../include/f_defs.h
vc_dv.o : ../include/f_defs.h
xk_wk_collect.o : ../include/f_defs.h
zstar_eu.o : ../include/f_defs.h
zstar_eu_us.o : ../include/f_defs.h

1
PP/make.depend
View File

@ -108,7 +108,6 @@ dipole.o : ../Modules/mp.o
dipole.o : ../Modules/mp_global.o
dipole.o : ../PW/para.o
dipole.o : ../PW/pwcom.o
dipole.o : ../PW/scf_mod.o
do_initial_state.o : ../Modules/cell_base.o
do_initial_state.o : ../Modules/io_global.o
do_initial_state.o : ../Modules/ions_base.o

4
PP/pw2casino.f90
View File

@ -69,7 +69,7 @@ SUBROUTINE compute_casino
nrxx, g, gg, ecutwfc, gcutm, nl, igtongl
USE klist , ONLY: nks, nelec, xk
USE lsda_mod, ONLY: lsda, nspin
USE scf, ONLY: rho, rho_core, rhog_core
USE scf, ONLY: rho, rho_core, rhog_core, tauk
USE vlocal, ONLY: vloc, vnew, strf
USE wvfct, ONLY: npw, npwx, nbnd, gamma_only, igk, g2kin, wg, et
USE uspp, ONLY: nkb, vkb, dvan
@ -219,7 +219,7 @@ SUBROUTINE compute_casino
!
! compute hartree and xc contribution
!
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, &
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, tauk, &
ehart, etxc, vtxc, etotefield, charge, vnew )
!
etot=(ek + (etxc-etxcc)+ehart+eloc+enl+ewld)

142
PW/electrons.f90
View File

@ -109,8 +109,8 @@ SUBROUTINE electrons()
! ... the charge density and of the HXC-potential
!
type (scf_type) :: rhoin ! used to store rho_in of current/next iteration
COMPLEX(DP), ALLOCATABLE :: taukgnew(:,:)
REAL(DP), ALLOCATABLE :: tauknew(:,:)
COMPLEX(DP), ALLOCATABLE :: taukgin(:,:)
REAL(DP), ALLOCATABLE :: taukin(:,:)
!
! ... external functions
!
@ -208,6 +208,10 @@ SUBROUTINE electrons()
becstep(:,:,:) = becsum(:,:,:)
END IF
call create_scf_type ( rhoin )
IF ( dft_is_meta() ) then
ALLOCATE( taukin( nrxx, nspin ) )
ALLOCATE( taukgin( ngm, nspin ) )
END IF
!
DO idum = 1, niter
!
@ -239,9 +243,11 @@ SUBROUTINE electrons()
deband_hwf = delta_e()
!
! save input current density in rhoin
rhoin%of_r(:,:) = rho%of_r(:,:)
rhoin%of_g(:,:) = rho%of_g(:,:)
!rhoin = rho
rhoin = rho
if ( dft_is_meta() ) then
taukin(:,:) = tauk(:,:)
taukgin(:,:) = taukg(:,:)
end if
!
scf_step: DO
!
@ -278,74 +284,10 @@ SUBROUTINE electrons()
! PAW : sum_band DOES NOT compute new one-center charges
CALL sum_band()
!
! ... bring output charge density (now stored in R-space in rho%of_r)
! ... also to G-space (rho%of_g) for mixing
!
DO is = 1, nspin
!
psic(:) = rho%of_r(:,is)
!
CALL cft3( psic, nr1, nr2, nr3, nrx1, nrx2, nrx3, -1 )
!
rho%of_g(:,is) = psic(nl(:))
!
IF ( okvan .AND. tqr ) THEN
!
! ... in case the augmentation charges are computed in real space
! ... we apply an FFT filter to the density in real space to
! ... remove features that are not compatible with the FFT grid.
!
psic(:) = ( 0.D0, 0.D0 )
!
psic(nl(:)) = rho%of_g(:,is)
!
IF ( gamma_only ) psic(nlm(:)) = CONJG( rho%of_g(:,is) )
!
CALL cft3( psic, nr1, nr2, nr3, nrx1, nrx2, nrx3, 1 )
!
rho%of_r(:,is) = psic(:)
!
END IF
!
END DO
!
! ... the same for tauk -> rhognew
!
IF ( dft_is_meta()) THEN
ALLOCATE( taukgnew( ngm, nspin ) )
DO is = 1, nspin
!
psic(:) = tauk(:,is)
!
CALL cft3( psic, nr1, nr2, nr3, nrx1, nrx2, nrx3, -1 )
!
taukgnew(:,is) = psic(nl(:))
!
IF ( okvan .AND. tqr ) THEN
!
! ... in case the augmentation terms are computed in real space
! ... we apply an FFT filter to the density in real space to
! ... remove features that are not compatible with the FFT grid
!
psic(:) = ( 0.D0, 0.D0 )
!
psic(nl(:)) = taukgnew(:,is)
!
IF ( gamma_only ) psic(nlm(:)) = CONJG( taukgnew(:,is) )
!
CALL cft3( psic, nr1, nr2, nr3, nrx1, nrx2, nrx3, 1 )
!
tauk(:,is) = psic(:)
!
END IF
!
END DO
END IF
!
! ... the Harris-Weinert-Foulkes energy is computed here using only
! ... quantities obtained from the input density
!
hwf_energy = eband + deband_hwf + ( etxc - etxcc ) + ewld + ehart + demet
hwf_energy = eband + deband_hwf + (etxc - etxcc) + ewld + ehart + demet
!
IF ( lda_plus_u ) THEN
!
@ -390,7 +332,7 @@ SUBROUTINE electrons()
becnew(:,:,:) = becsum(:,:,:)
END IF
!
CALL mix_rho( rho%of_g, rhoin%of_g, taukgnew, taukg, becnew, becstep, &
CALL mix_rho( rho%of_g, rhoin%of_g, taukg, taukgin, becnew, becstep, &
nsnew, ns, mixing_beta, dr2, tr2_min, iter, nmix, conv_elec )
!
! ... if convergence is achieved or if the self-consistency error
@ -398,19 +340,11 @@ SUBROUTINE electrons()
! ... (tr2_min), rhoin%of_g and rho%of_g are unchanged: rhoin%of_g
! ... contains the input density and rho%of_g contains the output
! ... density, both in G-space.
! ... in the other cases rhoin%of_g now contains mixed charge density
! ... (the new input density) in G-space.
! ... In the other cases rhoin%of_g now contains mixed charge density
! ... (the new input density) in G-space. In this case
! NB: In this later case mix_rho leaves out-of-sync rho(in)%of_r and
! rho(in)%of_g. This should be fixed as it is error-prone.
!
IF ( conv_elec ) THEN
!
! ... if convergence is achieved, rho%of_g and rho%of_r contain
! ... the same charge density, one in G-space, the other in R-space
!
IF ( dft_is_meta() ) taukg(:,:) = taukgnew(:,:)
!
END IF
!
IF ( dft_is_meta() ) DEALLOCATE( taukgnew )
IF ( okpaw ) DEALLOCATE (becnew)
!
IF ( first .and. nat > 0) THEN
@ -438,11 +372,11 @@ SUBROUTINE electrons()
not_converged_electrons : &
IF ( .NOT. conv_elec ) THEN
!
! ... bring the mixed charge density (rho%of_g) from G- to
! ... R-space (rhoin%of_r)
! ... synchronize R- and G- components of the mixed charge density
!
DO is = 1, nspin
!
! use psic as working array
psic(:) = ( 0.D0, 0.D0 )
!
psic(nl(:)) = rhoin%of_g(:,is)
@ -455,30 +389,28 @@ SUBROUTINE electrons()
!
END DO
!
! the same for the kinetic energy density (tauknew)
! the same for the kinetic energy density (taukin)
!
IF ( dft_is_meta() ) THEN
ALLOCATE( tauknew( nrxx, nspin ) )
DO is = 1, nspin
!
psic(:) = ( 0.D0, 0.D0 )
!
psic(nl(:)) = taukg(:,is)
psic(nl(:)) = taukgin(:,is)
!
IF ( gamma_only ) psic(nlm(:)) = CONJG( taukg(:,is) )
IF ( gamma_only ) psic(nlm(:)) = CONJG( taukgin(:,is) )
!
CALL cft3( psic, nr1, nr2, nr3, nrx1, nrx2, nrx3, 1 )
!
tauknew(:,is) = psic(:)
taukin(:,is) = psic(:)
!
END DO
!
END IF
!
! ... no convergence yet: calculate new potential from mixed
! ... charge density (i.e. the new estimate)
!
CALL v_of_rho( rhoin%of_r, rhoin%of_g, rho_core, rhog_core, &
CALL v_of_rho( rhoin%of_r, rhoin%of_g, rho_core, rhog_core, taukin, &
ehart, etxc, vtxc, etotefield, charge, vr )
!
! ... estimate correction needed to have variational energy:
@ -490,16 +422,13 @@ SUBROUTINE electrons()
!
descf = delta_escf()
!
! ... now copy the mixed charge density in R-space and G-space in rho
! ... now copy the mixed charge density in R- and G-space in rho
!
!rho%of_r(:,:) = rhoin%of_r(:,:)
!rho%of_g(:,:) = rhoin%of_g(:,:)
rho = rhoin
!
IF ( dft_is_meta() ) THEN
tauk(:,:) = tauknew(:,:)
DEALLOCATE( tauknew )
END IF
if ( dft_is_meta() ) then
tauk(:,:) = taukin(:,:)
taukg(:,:) = taukgin(:,:)
end if
!
! ... compute PAW corrections to descf
IF (okpaw) THEN
@ -537,7 +466,7 @@ SUBROUTINE electrons()
!
vnew(:,:) = vr(:,:)
!
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, &
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, tauk, &
ehart, etxc, vtxc, etotefield, charge, vr )
!
vnew(:,:) = vr(:,:) - vnew(:,:)
@ -673,7 +602,7 @@ SUBROUTINE electrons()
#endif
correction1c = (deband_PAW + descf_PAW + e_PAW)
PRINT '(5x,A,f10.6,A)', 'PAW correction: ',correction1c, ', composed of: '
PRINT '(5x,A,f10.6,A,f10.6,A,f12.6)',&
PRINT '(5x,A,f10.6,A,f10.6,A,f10.6)',&
'de_band: ',deband_PAW,', de_scf: ',descf_PAW,', 1-center E: ', e_PAW
IF (really_do_paw) THEN
etot = etot + correction1c
@ -694,7 +623,7 @@ SUBROUTINE electrons()
IF ( first ) THEN
!
fock0 = exxenergy2()
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, &
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, tauk, &
ehart, etxc, vtxc, etotefield, charge, vr )
!
CALL set_vrs( vrs, vltot, vr, nrxx, nspin, doublegrid )
@ -824,6 +753,11 @@ SUBROUTINE electrons()
DEALLOCATE(deband_1ae_na, deband_1ps_na, descf_1ae_na, descf_1ps_na)
#endif
END IF
!
IF ( dft_is_meta() ) then
DEALLOCATE( taukgin )
DEALLOCATE( taukin )
ENDIF
call destroy_scf_type ( rhoin )
!
RETURN
@ -1049,7 +983,7 @@ SUBROUTINE electrons()
IF ( dft_is_meta() ) THEN
DO ipol = 1, nspin
delta_escf = delta_escf - &
SUM( (tauknew(:,ipol)-tauk(:,ipol) )*kedtaur(:,ipol))
SUM( (taukin(:,ipol)-tauk(:,ipol) )*kedtaur(:,ipol))
END DO
END IF
!

6
PW/make.depend
View File

@ -814,7 +814,9 @@ paw_xc.o : exx.o
paw_xc.o : noncol.o
paw_xc.o : pwcom.o
potinit.o : ../Modules/cell_base.o
potinit.o : ../Modules/constants.o
potinit.o : ../Modules/control_flags.o
potinit.o : ../Modules/functionals.o
potinit.o : ../Modules/grid_paw_variables.o
potinit.o : ../Modules/io_files.o
potinit.o : ../Modules/io_global.o
@ -948,6 +950,7 @@ read_conf_from_file.o : pw_restart.o
read_conf_from_file.o : pwcom.o
read_file.o : ../Modules/cell_base.o
read_file.o : ../Modules/constants.o
read_file.o : ../Modules/grid_paw_variables.o
read_file.o : ../Modules/io_files.o
read_file.o : ../Modules/ions_base.o
read_file.o : ../Modules/kind.o
@ -956,9 +959,11 @@ read_file.o : ../Modules/uspp.o
read_file.o : ../Modules/wavefunctions.o
read_file.o : ../Modules/xml_io_base.o
read_file.o : buffers.o
read_file.o : grid_paw_routines.o
read_file.o : noncol.o
read_file.o : pw_restart.o
read_file.o : pwcom.o
read_file.o : rad_paw_routines.o
read_file.o : scf_mod.o
read_ncpp.o : ../Modules/functionals.o
read_ncpp.o : ../Modules/kind.o
@ -1243,6 +1248,7 @@ sum_band.o : ../Modules/wavefunctions.o
sum_band.o : buffers.o
sum_band.o : noncol.o
sum_band.o : pwcom.o
sum_band.o : realus.o
sum_band.o : scf_mod.o
sumkg.o : ../Modules/kind.o
sumkt.o : ../Modules/kind.o

19
PW/potinit.f90
View File

@ -24,6 +24,7 @@ SUBROUTINE potinit()
! ... In all cases the scf potential is recalculated and saved in vr
!
USE kinds, ONLY : DP
USE constants, ONLY : pi
USE io_global, ONLY : stdout
USE cell_base, ONLY : alat, omega
USE ions_base, ONLY : nat, ityp, ntyp => nsp
@ -34,8 +35,9 @@ SUBROUTINE potinit()
ngm, gstart, nl, g, gg
USE gsmooth, ONLY : doublegrid
USE control_flags, ONLY : lscf
USE scf, ONLY : rho, rho_core, rhog_core, &
USE scf, ONLY : rho, rho_core, rhog_core, tauk, taukg, &
vltot, vr, vrs
USE funct, ONLY : dft_is_meta
USE wavefunctions_module, ONLY : psic
USE ener, ONLY : ehart, etxc, vtxc
USE ldaU, ONLY : niter_with_fixed_ns
@ -63,6 +65,7 @@ SUBROUTINE potinit()
!
REAL(DP) :: charge ! the starting charge
REAL(DP) :: etotefield !
REAL(DP) :: fact
INTEGER :: is, ios
INTEGER :: ldim ! integer variable for I/O control
LOGICAL :: exst
@ -217,9 +220,21 @@ SUBROUTINE potinit()
!
END DO
!
if ( dft_is_meta()) then
! ... define a starting (TF) guess for tauk and taukg
fact = (3.d0*pi*pi)**(2.0/3.0)
DO is = 1, nspin
tauk(:,is) = fact * abs(rho%of_r(:,is)*nspin)**(5.0/3.0)/nspin
psic(:) = tauk(:,is)
CALL cft3( psic, nr1, nr2, nr3, nrx1, nrx2, nrx3, -1 )
taukg(:,is) = psic(nl(:))
END DO
!
end if
!
! ... compute the potential and store it in vr
!
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, &
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, tauk, &
ehart, etxc, vtxc, etotefield, charge, vr )
!
! ... define the total local potential (external+scf)

4
PW/read_file.f90
View File

@ -31,7 +31,7 @@ SUBROUTINE read_file()
nl, gstart
USE gsmooth, ONLY : ngms, nls, nrx1s, nr1s, nr2s, nr3s
USE spin_orb, ONLY : so, lspinorb
USE scf, ONLY : rho, rho_core, rhog_core, vr
USE scf, ONLY : rho, rho_core, rhog_core, tauk, vr
USE wavefunctions_module, ONLY : psic
USE vlocal, ONLY : strf
USE io_files, ONLY : tmp_dir, prefix, iunpun, nwordwfc, iunwfc
@ -212,7 +212,7 @@ SUBROUTINE read_file()
!
! ... recalculate the potential
!
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, &
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, tauk, &
ehart, etxc, vtxc, etotefield, charge, vr )
!
! ... reads the wavefunctions and writes them in 'distributed' form

4
PW/restart_in_ions.f90
View File

@ -18,7 +18,7 @@ subroutine restart_in_ions (iter, ik_, dr2)
nrxx
USE klist, ONLY: nks
USE lsda_mod, ONLY: nspin
USE scf, ONLY : rho, rho_core, rhog_core
USE scf, ONLY : rho, rho_core, rhog_core, tauk
USE control_flags, ONLY: restart, tr2, ethr
USE vlocal, ONLY: vnew
USE wvfct, ONLY: nbnd, et
@ -80,7 +80,7 @@ subroutine restart_in_ions (iter, ik_, dr2)
! recalculate etxc, vtxc, ehart, needed by stress calculation
!
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, &
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, tauk, &
ehart, etxc, vtxc, etotefield, charge, psic )
!
! restart procedure completed

58
PW/sum_band.f90
View File

@ -21,13 +21,15 @@ SUBROUTINE sum_band()
USE control_flags, ONLY : diago_full_acc
USE cell_base, ONLY : at, bg, omega, tpiba
USE ions_base, ONLY : nat, ntyp => nsp, ityp
USE gvect, ONLY : nr1, nr2, nr3, nrx1, nrx2, nrx3, nrxx, ngm, g
USE gvect, ONLY : nr1, nr2, nr3, nrx1, nrx2, nrx3, nrxx, &
ngm, g, nl, nlm
USE gsmooth, ONLY : nls, nlsm, nr1s, nr2s, nr3s, &
nrx1s, nrx2s, nrx3s, nrxxs, doublegrid
USE klist, ONLY : nks, nkstot, wk, xk, ngk
USE ldaU, ONLY : lda_plus_U
USE lsda_mod, ONLY : lsda, nspin, current_spin, isk
USE scf, ONLY : rho, tauk
USE realus, ONLY : tqr
USE scf, ONLY : rho, tauk, taukg
USE symme, ONLY : nsym, s, ftau
USE io_files, ONLY : iunwfc, nwordwfc, iunigk
USE buffers, ONLY : get_buffer
@ -185,7 +187,57 @@ SUBROUTINE sum_band()
!
#endif
if (dft_is_meta() ) deallocate (kplusg)
!
!
! ... synchronize rho%of_g to the calculated rho%of_r
!
DO is = 1, nspin
!
! use psic as work array
psic(:) = rho%of_r(:,is)
CALL cft3( psic, nr1, nr2, nr3, nrx1, nrx2, nrx3, -1 )
rho%of_g(:,is) = psic(nl(:))
!
IF ( okvan .AND. tqr ) THEN
! ... in case the augmentation charges are computed in real space
! ... we apply an FFT filter to the density in real space to
! ... remove features that are not compatible with the FFT grid.
!
psic(:) = ( 0.D0, 0.D0 )
psic(nl(:)) = rho%of_g(:,is)
IF ( gamma_only ) psic(nlm(:)) = CONJG( rho%of_g(:,is) )
CALL cft3( psic, nr1, nr2, nr3, nrx1, nrx2, nrx3, 1 )
rho%of_r(:,is) = psic(:)
!
END IF
!
END DO
!
! ... the same for tauk and taukg
!
IF ( dft_is_meta()) THEN
DO is = 1, nspin
!
! use psic as work array
psic(:) = tauk(:,is)
CALL cft3( psic, nr1, nr2, nr3, nrx1, nrx2, nrx3, -1 )
taukg(:,is) = psic(nl(:))
!
IF ( okvan .AND. tqr ) THEN
! ... in case the augmentation charges are computed in real space
! ... we apply an FFT filter to the density in real space to
! ... remove features that are not compatible with the FFT grid.
!
psic(:) = ( 0.D0, 0.D0 )
psic(nl(:)) = taukg(:,is)
IF ( gamma_only ) psic(nlm(:)) = CONJG( taukg(:,is) )
CALL cft3( psic, nr1, nr2, nr3, nrx1, nrx2, nrx3, 1 )
tauk(:,is) = psic(:)
!
END IF
!
END DO
END IF
!
CALL stop_clock( 'sum_band' )
!
RETURN

4
PW/update_pot.f90
View File

@ -192,7 +192,7 @@ SUBROUTINE extrapolate_charge( rho_extr )
USE gvect, ONLY : nrxx, ngm, g, gg, gstart, nr1, nr2, nr3, &
nl, eigts1, eigts2, eigts3, nrx1, nrx2, nrx3
USE lsda_mod, ONLY : lsda, nspin
USE scf, ONLY : rho, rho_core, rhog_core, vr
USE scf, ONLY : rho, rho_core, rhog_core, tauk, vr
USE wavefunctions_module, ONLY : psic
USE control_flags, ONLY : alpha0, beta0
USE ener, ONLY : ehart, etxc, vtxc
@ -378,7 +378,7 @@ SUBROUTINE extrapolate_charge( rho_extr )
!
END DO
!
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, &
CALL v_of_rho( rho%of_r, rho%of_g, rho_core, rhog_core, tauk, &
ehart, etxc, vtxc, etotefield, charge, vr )
!
IF ( ABS( charge - nelec ) / charge > 1.D-7 ) THEN

18
PW/v_of_rho.f90
View File

@ -6,7 +6,7 @@
! or http://www.gnu.org/copyleft/gpl.txt .
!
!----------------------------------------------------------------------------
SUBROUTINE v_of_rho( rho, rhog, rho_core, rhog_core, &
SUBROUTINE v_of_rho( rho, rhog, rho_core, rhog_core, tauk, &
ehart, etxc, vtxc, etotefield, charge, v )
!----------------------------------------------------------------------------
!
@ -23,7 +23,7 @@ SUBROUTINE v_of_rho( rho, rhog, rho_core, rhog_core, &
!
IMPLICIT NONE
!
REAL(DP), INTENT(IN) :: rho(nrxx,nspin), rho_core(nrxx)
REAL(DP), INTENT(IN) :: rho(nrxx,nspin), rho_core(nrxx), tauk(nrxx,nspin)
! input: the valence charge
! input: the core charge
COMPLEX(DP), INTENT(IN) :: rhog(ngm,nspin), rhog_core(ngm)
@ -43,7 +43,7 @@ SUBROUTINE v_of_rho( rho, rhog, rho_core, rhog_core, &
! ... calculate exchange-correlation potential
!
if (dft_is_meta()) then
call v_xc_meta( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
call v_xc_meta( rho, rhog, rho_core, rhog_core, tauk, etxc, vtxc, v )
else
CALL v_xc( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
endif
@ -70,7 +70,7 @@ SUBROUTINE v_of_rho( rho, rhog, rho_core, rhog_core, &
!
END SUBROUTINE v_of_rho
!
SUBROUTINE v_xc_meta( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
SUBROUTINE v_xc_meta( rho, rhog, rho_core, rhog_core, tauk, etxc, vtxc, v )
!----------------------------------------------------------------------------
!
! ... Exchange-Correlation potential Vxc(r) from n(r)
@ -86,7 +86,7 @@ SUBROUTINE v_xc_meta( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
!
IMPLICIT NONE
!
REAL(DP), INTENT(IN) :: rho(nrxx,nspin), rho_core(nrxx)
REAL(DP), INTENT(IN) :: rho(nrxx,nspin), rho_core(nrxx), tauk(nrxx,nspin)
! the valence charge
! the core charge
COMPLEX(DP), INTENT(IN) :: rhog(ngm,nspin), rhog_core(ngm)
@ -123,7 +123,7 @@ SUBROUTINE v_xc_meta( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
rhoneg = 0.D0
!
! IF (get_igcx()==7.AND.get_igcc()==6) THEN
call v_xc_tpss( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
call v_xc_tpss( rho, rhog, rho_core, rhog_core, tauk, etxc, vtxc, v )
! ELSE
! CALL errore('v_xc_meta','wrong igcx and/or igcc',1)
! ENDIF
@ -131,7 +131,7 @@ SUBROUTINE v_xc_meta( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
RETURN
END SUBROUTINE v_xc_meta
!
SUBROUTINE v_xc_tpss( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
SUBROUTINE v_xc_tpss( rho, rhog, rho_core, rhog_core, tauk, etxc, vtxc, v )
! ===================
!--------------------------------------------------------------------
! use gvecp, only: ng => ngm
@ -139,14 +139,14 @@ SUBROUTINE v_xc_tpss( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
USE gvect, ONLY : nrxx, nrx1,nrx2,nrx3,nr1,nr2,nr3, &
g,nl,ngm
USE gsmooth, ONLY : nr1s, nr2s, nr3s, nrx1s, nrx2s, nrx3s, nrxxs
USE scf, ONLY : kedtau,tauk, kedtaur
USE scf, ONLY : kedtau, kedtaur
USE lsda_mod, ONLY : nspin
USE cell_base, ONLY : omega, alat
USE constants, ONLY : e2
IMPLICIT NONE
!
! input
REAL(DP),INTENT(IN) :: rho(nrxx,nspin), rho_core(nrxx)
REAL(DP),INTENT(IN) :: rho(nrxx,nspin), rho_core(nrxx), tauk(nrxx,nspin)
COMPLEX(DP),INTENT(IN) :: rhog(ngm,nspin), rhog_core(ngm)
REAL(DP),INTENT(OUT) :: etxc, vtxc, v(nrxx,nspin)
! integer nspin , nnr