TABLE OF CONTENTS
ABINIT/newfermie1 [ Functions ]
NAME
newfermie1
FUNCTION
This routine computes the derivative of the fermi energy wrt the active perturbation for use in evaluating the edocc term and active subspace contribution to the first-order wavefunctions in the case of metals. This is presently used only for the strain and magnetic field perturbations, and only for Q = 0.
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (DRH) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
cplex: if 1, real space 1-order functions on FFT grid are REAL, if 2, COMPLEX fe1fixed=fixed contribution to the first-order Fermi energy ipert=index of perturbation istep=index of the number of steps in the routine scfcv ixc= choice of exchange-correlation scheme mpi_atmtab(:)=--optional-- indexes of the atoms treated by current proc comm_atom=--optional-- MPI communicator over atoms my_natom=number of atoms treated by current processor natom=number of atoms nfft=(effective) number of FFT grid points (for this processor) nfftot= total number of FFT grid points nhatfermi(nfft,nspden)=fermi-level compensation charge density (PAW only) nspden=number of spin-density components ntypat=number of atom types occopt=option for occupancies paw_an(natom) <type(paw_an_type)>=paw arrays for 0th-order quantities given on angular mesh paw_an1(natom) <type(paw_an_type)>=paw arrays for 1st-order quantities given on angular mesh paw_ij1(natom) <type(paw_ij_type)>=(1st-order) paw arrays given on (i,j) channels pawang <type(pawang_type)>=paw angular mesh and related data pawnzlm=-- PAW only -- option for the computation of non-zero lm moments of the on-sites densities pawrad(ntypat) <type(pawrad_type)>=paw radial mesh and related data pawrhoij1(natom) <type(pawrhoij_type)>= paw rhoij 1st-order occupancies pawrhoijfermi(natom) <type(pawrhoij_type)>=paw rhoij occupancies at Fermi level pawtab(ntypat) <type(pawtab_type)>=paw tabulated starting data pawxcdev=Choice of XC development (0=no dev. (use of angular mesh) ; 1 or 2=dev. on moments) xclevel= XC functional level prtvol=control print volume and debugging output rhorfermi(nfft,nspden)=fermi-level electronic density ucvol=unit cell volume in bohr**3 usepaw=1 if PAW is activated usexcnhat= -PAW only- flag controling use of compensation density in Vxc vtrial1(cplex*nfft,nspden)=1-st order potential vxc1(cplex*nfft,nspden)=1-st order XC potential
OUTPUT
(see side effects)
SIDE EFFECTS
fermie1=derivative of fermi energy wrt perturbation at input : old value at output : updated value
PARENTS
dfpt_scfcv
CHILDREN
dotprod_vn,free_my_atmtab,get_my_atmtab,pawdfptenergy,wrtout
SOURCE
74 #if defined HAVE_CONFIG_H 75 #include "config.h" 76 #endif 77 78 #include "abi_common.h" 79 80 81 subroutine newfermie1(cplex,fermie1,fe1fixed,ipert,istep,ixc,my_natom,natom,nfft,nfftot,& 82 & nhatfermi,nspden,ntypat,occopt,paw_an,paw_an1,paw_ij1,pawang,pawnzlm,pawrad,& 83 & pawrhoij1,pawrhoijfermi,pawtab,pawxcdev,prtvol,rhorfermi,& 84 & ucvol,usepaw,usexcnhat,vtrial1,vxc1,xclevel,& 85 & mpi_atmtab,comm_atom) ! optional arguments (parallelism) 86 87 use defs_basis 88 use m_profiling_abi 89 use m_errors 90 use m_xmpi 91 92 use m_pawang, only : pawang_type 93 use m_pawrad, only : pawrad_type 94 use m_pawtab, only : pawtab_type 95 use m_paw_an, only : paw_an_type 96 use m_paw_ij, only : paw_ij_type 97 use m_pawrhoij, only : pawrhoij_type 98 use m_paral_atom, only : get_my_atmtab, free_my_atmtab 99 use m_cgtools, only : dotprod_vn 100 101 !This section has been created automatically by the script Abilint (TD). 102 !Do not modify the following lines by hand. 103 #undef ABI_FUNC 104 #define ABI_FUNC 'newfermie1' 105 use interfaces_14_hidewrite 106 use interfaces_65_paw 107 !End of the abilint section 108 109 implicit none 110 111 !Arguments ------------------------------- 112 !scalars 113 integer,intent(in) :: cplex,ipert,istep,ixc,my_natom,natom,nfft,nfftot,nspden,ntypat 114 integer,intent(in) :: occopt,pawnzlm,pawxcdev,prtvol,usepaw,usexcnhat,xclevel 115 integer,optional,intent(in) :: comm_atom 116 real(dp),intent(in) :: fe1fixed,ucvol 117 real(dp),intent(inout) :: fermie1 118 type(pawang_type),intent(in) :: pawang 119 !arrays 120 integer,optional,target,intent(in) :: mpi_atmtab(:) 121 real(dp),intent(in) :: rhorfermi(nfft,nspden),vtrial1(cplex*nfft,nspden) 122 real(dp),intent(in) :: nhatfermi(:,:),vxc1(:,:) 123 type(paw_an_type),intent(in) :: paw_an(my_natom*usepaw) 124 type(paw_an_type),intent(inout) :: paw_an1(my_natom*usepaw) 125 type(paw_ij_type),intent(inout) :: paw_ij1(my_natom*usepaw) 126 type(pawrad_type),intent(in) :: pawrad(ntypat*usepaw) 127 type(pawrhoij_type),intent(in) :: pawrhoij1(my_natom*usepaw),pawrhoijfermi(my_natom*usepaw) 128 type(pawtab_type),intent(in) :: pawtab(ntypat*usepaw) 129 130 !Local variables------------------------------- 131 !scalars 132 integer :: ipert0,my_comm_atom,nzlmopt,nzlmopt_fermi,option,pawprtvol 133 logical :: my_atmtab_allocated,paral_atom 134 real(dp) :: doti,fe1_scf,fe1_tmp,fermie1_new,fermie1rs 135 character(len=500) :: msg 136 !arrays 137 integer, pointer :: my_atmtab(:) 138 real(dp) :: fe1_paw(2) 139 real(dp), allocatable :: rhor_nonhat(:,:),vtrial1_novxc(:,:) 140 141 ! ********************************************************************* 142 143 !Tests 144 if (cplex==2) then 145 msg='Not compatible with cplex=2!' 146 MSG_BUG(msg) 147 end if 148 if (usepaw==1.and.usexcnhat==0.and.(size(nhatfermi)<=0.or.size(vxc1)<=0)) then 149 msg='Should have nhatfermi and vxc1 allocated with usexcnhat=0!' 150 MSG_BUG(msg) 151 end if 152 153 !Set up parallelism over atoms 154 paral_atom=(present(comm_atom).and.(my_natom/=natom)) 155 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 156 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 157 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 158 159 if(occopt>=3 .and. occopt <=8) then 160 161 ! The product of the current trial potential and the so-called Fermi level 162 ! density is integrated to give the local potential contributions to the 163 ! first-order Fermi level. 164 option=1 165 if (usepaw==1.and.usexcnhat==0) then 166 ABI_ALLOCATE(rhor_nonhat,(nfft,nspden)) 167 ABI_ALLOCATE(vtrial1_novxc,(nfft,nspden)) 168 rhor_nonhat(1:nfft,1:nspden)=rhorfermi(1:nfft,1:nspden)-nhatfermi(1:nfft,1:nspden) 169 vtrial1_novxc(1:nfft,1:nspden)=vtrial1(1:nfft,1:nspden)-vxc1(1:nfft,1:nspden) 170 call dotprod_vn(cplex,rhor_nonhat,fe1_scf,doti,nfft,nfftot,& 171 & nspden,option,vtrial1,ucvol) 172 call dotprod_vn(cplex,nhatfermi,fe1_tmp,doti,nfft,nfftot,& 173 & nspden,option,vtrial1_novxc,ucvol) 174 fe1_scf=fe1_scf+fe1_tmp 175 ABI_DEALLOCATE(rhor_nonhat) 176 ABI_DEALLOCATE(vtrial1_novxc) 177 else 178 call dotprod_vn(cplex,rhorfermi,fe1_scf,doti,nfft,nfftot,& 179 & nspden,option,vtrial1,ucvol) 180 end if 181 182 fe1_paw(:)=zero 183 ! PAW on-site contribution (use Fermi level occupation matrix) 184 if (usepaw==1) then 185 ipert0=0;pawprtvol=0 186 nzlmopt=0;if (istep>1) nzlmopt=pawnzlm 187 if (istep==1.and.pawnzlm>0) nzlmopt=-1 188 nzlmopt_fermi=0;if (pawnzlm>0) nzlmopt_fermi=-1 189 call pawdfptenergy(fe1_paw,ipert,ipert0,ixc,my_natom,natom,ntypat,nzlmopt,& 190 & nzlmopt_fermi,paw_an,paw_an1,paw_ij1,pawang,pawprtvol,pawrad,& 191 & pawrhoij1,pawrhoijfermi,pawtab,pawxcdev,xclevel,& 192 & mpi_atmtab=my_atmtab, comm_atom=my_comm_atom) 193 end if 194 195 ! The fixed contributions consisting of non-local potential and kinetic terms 196 ! are added 197 fermie1_new=fe1fixed+fe1_scf+fe1_paw(1) 198 fermie1rs=(fermie1-fermie1_new)**2 199 fermie1=fermie1_new 200 201 if(prtvol>=10)then 202 write(msg, '(a,i5,2es18.8)' ) ' fermie1, residual squared',istep,fermie1,fermie1rs 203 call wrtout(std_out,msg,'COLL') 204 end if 205 206 else 207 fermie1=zero 208 end if 209 210 !Destroy atom table used for parallelism 211 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 212 213 end subroutine newfermie1