TABLE OF CONTENTS
ABINIT/dfpt_dyfro [ Functions ]
NAME
dfpt_dyfro
FUNCTION
Compute the different parts of the frozen-wavefunction part of the dynamical matrix, except the non-local one, computed previously. Also (when installed) symmetrize the different part and their sum.
COPYRIGHT
Copyright (C) 2000-2018 ABINIT group (XG,MT) 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
atindx1(natom)=index table for atoms, inverse of atindx dyfr_cplex=1 if dyfrnl is real, 2 if it is complex dyfr_nondiag=1 if dyfrnl and dyfrwf are non diagonal with respect to atoms; 0 otherwise gmet(3,3)=reciprocal space metric (bohr^-2) gprimd(3,3)=dimensional primitive translations for reciprocal space (bohr**-1) gsqcut=cutoff on G^2 based on ecut indsym(4,nsym,natom)=index showing transformation of atom labels under symmetry operations (computed in symatm) mgfft=maximum size of 1D FFTs mpi_enreg=information about MPI parallelization mqgrid=dimensioned number of q grid points for local psp spline natom=number of atoms in unit cell nattyp(ntypat)=number of atoms of each type nfft=(effective) number of FFT grid points (for this processor) ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft nspden=number of spin-density components nsym=number of symmetries in space group ntypat=number of types of atoms n1xccc=dimension of xccc1d ; 0 if no XC core correction is used n3xccc=dimension of the xccc3d array (0 or nfft). psps <type(pseudopotential_type)>=variables related to pseudopotentials pawtab(ntypat*usepaw) <type(pawtab_type)>=paw tabulated starting data ph1d(2,3*(2*mgfft+1)*natom)=one-dimensional structure factor information qgrid(mqgrid)=q point array for local psp spline fits qphon(3)=wavevector of the phonon rhog(2,nfft)=electron density in G space rprimd(3,3)=dimensional primitive translation vectors (bohr) symq(4,2,nsym)=1 if symmetry preserves present qpoint. From littlegroup_q symrec(3,3,nsym)=symmetries in reciprocal space typat(natom)=integer type for each atom in cell ucvol=unit cell volume (bohr**3). usepaw= 0 for non paw calculation; =1 for paw calculation vlspl(mqgrid,2,ntypat)=q^2 v(q) spline for each type of atom. vxc(nfft,nspden)=exchange-correlation potential (hartree) in real space--only used when n1xccc/=0 xcccrc(ntypat)=XC core correction cutoff radius (bohr) for each atom type xccc1d(n1xccc,6,ntypat)=1D core charge function and five derivatives, for each type of atom, from psp xccc3d(n3xccc)=3D core electron density for XC core correction, bohr^-3 xred(3,natom)=reduced coordinates for atoms in unit cell
OUTPUT
dyfrlo(3,3,natom)=frozen wavefunctions part of the dynamical matrix (local only) dyfrwf(dyfr_cplex,3,3,natom,1+(natom-1)*dyfr_nondiag)= frozen wavefunctions part of the dynamical matrix (local + non-local) If NCPP, it depends on one atom If PAW, it depends on two atoms dyfrxc(3,3,natom)=frozen wavefunctions part of the dynamical matrix (non-linear xc core correction)
SIDE EFFECTS
Input/Output dyfrnl(dyfr_cplex,3,3,natom,1+(natom-1)*dyfr_nondiag)= frozen wavefunctions part of the dynamical matrix (non-local only) If NCPP, it depends on one atom If PAW, it depends on two atoms
PARENTS
respfn
CHILDREN
atm2fft,dfpt_sydy,fourdp,mkcore,mklocl_recipspace,timab,zerosym
SOURCE
90 #if defined HAVE_CONFIG_H 91 #include "config.h" 92 #endif 93 94 #include "abi_common.h" 95 96 subroutine dfpt_dyfro(atindx1,dyfrnl,dyfrlo,dyfrwf,dyfrxc,dyfr_cplex,dyfr_nondiag,& 97 & gmet,gprimd,gsqcut,indsym,mgfft,mpi_enreg,mqgrid,natom,nattyp,& 98 & nfft,ngfft,nspden,nsym,ntypat,n1xccc,n3xccc,paral_kgb,psps,pawtab,ph1d,qgrid,& 99 & qphon,rhog,rprimd,symq,symrec,typat,ucvol,usepaw,vlspl,vxc,& 100 & xcccrc,xccc1d,xccc3d,xred) 101 102 use defs_basis 103 use defs_abitypes 104 use m_profiling_abi 105 use m_errors 106 107 use defs_datatypes, only : pseudopotential_type 108 use m_pawtab, only : pawtab_type 109 use m_fft, only : zerosym 110 111 !This section has been created automatically by the script Abilint (TD). 112 !Do not modify the following lines by hand. 113 #undef ABI_FUNC 114 #define ABI_FUNC 'dfpt_dyfro' 115 use interfaces_18_timing 116 use interfaces_53_ffts 117 use interfaces_56_xc 118 use interfaces_64_psp 119 use interfaces_67_common 120 use interfaces_72_response, except_this_one => dfpt_dyfro 121 !End of the abilint section 122 123 implicit none 124 125 !Arguments ------------------------------------ 126 !scalars 127 integer,intent(in) :: dyfr_cplex,dyfr_nondiag,mgfft,mqgrid,n1xccc,n3xccc,natom,nfft,nspden 128 integer,intent(in) :: nsym,ntypat,paral_kgb,usepaw 129 real(dp),intent(in) :: gsqcut,ucvol 130 type(pseudopotential_type),intent(in) :: psps 131 type(MPI_type),intent(in) :: mpi_enreg 132 !arrays 133 integer,intent(in) :: atindx1(natom),indsym(4,nsym,natom),nattyp(ntypat) 134 integer,intent(in) :: ngfft(18),symq(4,2,nsym),symrec(3,3,nsym),typat(natom) 135 real(dp),intent(in) :: gmet(3,3),gprimd(3,3),ph1d(2,3*(2*mgfft+1)*natom) 136 real(dp),intent(in) :: qgrid(mqgrid),qphon(3),rhog(2,nfft),rprimd(3,3) 137 real(dp),intent(in) :: vlspl(mqgrid,2,ntypat),vxc(nfft,nspden) 138 real(dp),intent(in) :: xccc1d(n1xccc,6,ntypat),xcccrc(ntypat),xred(3,natom) 139 real(dp),intent(inout) :: dyfrnl(dyfr_cplex,3,3,natom,1+(natom-1)*dyfr_nondiag),xccc3d(n3xccc) 140 real(dp),intent(out) :: dyfrlo(3,3,natom),dyfrwf(dyfr_cplex,3,3,natom,1+(natom-1)*dyfr_nondiag) 141 real(dp),intent(inout) :: dyfrxc(3,3,natom) !vz_i 142 type(pawtab_type),intent(in) :: pawtab(ntypat*usepaw) 143 144 !Local variables------------------------------- 145 !scalars 146 logical, parameter :: do_final_sym=.true. 147 integer :: iatom,jatom,n1,n2,n3,optatm,optdyfr,opteltfr,optgr,option 148 integer :: optn,optn2,optstr,optv 149 real(dp) :: eei 150 !arrays 151 integer :: qprtrb(3) 152 real(dp) :: dummy6(6),dum_strn(6),dum_strv(6) 153 real(dp) :: tsec(2),vprtrb(2) 154 real(dp) :: dum_atmrho(0),dum_atmvloc(0),dum_gauss(0),dum_grn(0),dum_grv(0),dum_eltfrxc(0) 155 real(dp),allocatable :: dyfrlo_tmp1(:,:,:),dyfrlo_tmp2(:,:,:,:,:),dyfrsym_tmp(:,:,:,:,:) 156 real(dp),allocatable :: gr_dum(:,:),v_dum(:),vxctotg(:,:) 157 158 ! ************************************************************************* 159 160 if(nspden==4)then 161 MSG_WARNING('dfpt_dyfro : DFPT with nspden=4 works at the moment just for insulators and norm-conserving psp!') 162 end if 163 164 n1=ngfft(1); n2=ngfft(2); n3=ngfft(3) 165 166 if (usepaw==1 .or. psps%nc_xccc_gspace==1) then 167 168 ! PAW or NC with nc_xccc_gspace: compute local psp and core charge contribs together 169 ! in reciprocal space 170 ! ----------------------------------------------------------------------- 171 call timab(563,1,tsec) 172 if (n3xccc>0) then 173 ABI_ALLOCATE(v_dum,(nfft)) 174 ABI_ALLOCATE(vxctotg,(2,nfft)) 175 v_dum(:)=vxc(:,1);if (nspden>=2) v_dum(:)=0.5_dp*(v_dum(:)+vxc(:,2)) 176 call fourdp(1,vxctotg,v_dum,-1,mpi_enreg,nfft,ngfft,paral_kgb,0) 177 call zerosym(vxctotg,2,ngfft(1),ngfft(2),ngfft(3),& 178 & comm_fft=mpi_enreg%comm_fft,distribfft=mpi_enreg%distribfft) 179 ABI_DEALLOCATE(v_dum) 180 else 181 ABI_ALLOCATE(vxctotg,(0,0)) 182 end if 183 optatm=0;optdyfr=1;optgr=0;optstr=0;optv=1;optn=n3xccc/nfft;optn2=1;opteltfr=0 184 call atm2fft(atindx1,dum_atmrho,dum_atmvloc,dyfrxc,dyfrlo,dum_eltfrxc,& 185 & dum_gauss,gmet,gprimd,dum_grn,dum_grv,gsqcut,mgfft,mqgrid,natom,nattyp,nfft,ngfft,ntypat,& 186 & optatm,optdyfr,opteltfr,optgr,optn,optn2,optstr,optv,psps,pawtab,ph1d,qgrid,qprtrb,rhog,& 187 & dum_strn,dum_strv,ucvol,usepaw,vxctotg,vxctotg,vxctotg,vprtrb,vlspl) 188 ABI_DEALLOCATE(vxctotg) 189 if (n3xccc==0) dyfrxc=zero 190 else 191 192 ! Norm-conserving: compute local psp contribution in reciprocal space 193 ! and core charge contribution in real space 194 ! ----------------------------------------------------------------------- 195 option=4 196 ABI_ALLOCATE(dyfrlo_tmp1,(3,3,natom)) 197 ABI_ALLOCATE(gr_dum,(3,natom)) 198 ABI_ALLOCATE(v_dum,(nfft)) 199 call mklocl_recipspace(dyfrlo_tmp1,eei,gmet,gprimd,& 200 & gr_dum,gsqcut,dummy6,mgfft,mpi_enreg,mqgrid,natom,nattyp,nfft,ngfft,& 201 & ntypat,option,paral_kgb,ph1d,qgrid,qprtrb,rhog,ucvol,vlspl,vprtrb,v_dum) 202 do iatom=1,natom 203 ! Reestablish correct order of atoms 204 dyfrlo(1:3,1:3,atindx1(iatom))=dyfrlo_tmp1(1:3,1:3,iatom) 205 end do 206 ABI_DEALLOCATE(dyfrlo_tmp1) 207 ABI_DEALLOCATE(v_dum) 208 if(n1xccc/=0)then 209 call mkcore(dummy6,dyfrxc,gr_dum,mpi_enreg,natom,nfft,nspden,ntypat,& 210 & n1,n1xccc,n2,n3,option,rprimd,typat,ucvol,vxc,xcccrc,xccc1d,xccc3d,xred) 211 end if 212 ABI_DEALLOCATE(gr_dum) 213 end if 214 215 !Symmetrize dynamical matrix explicitly for given space group: 216 217 !Symmetrize local part of the dynamical matrix dyfrlo: 218 ABI_ALLOCATE(dyfrsym_tmp,(1,3,3,natom,1)) 219 ABI_ALLOCATE(dyfrlo_tmp2,(1,3,3,natom,1)) 220 dyfrsym_tmp(1,:,:,:,1)=dyfrlo(:,:,:) 221 call dfpt_sydy(1,dyfrsym_tmp,indsym,natom,0,nsym,qphon,dyfrlo_tmp2,symq,symrec) 222 dyfrlo(:,:,:)=dyfrlo_tmp2(1,:,:,:,1) 223 if (do_final_sym) then 224 dyfrsym_tmp(1,:,:,:,1)=dyfrxc(:,:,:) 225 call dfpt_sydy(1,dyfrsym_tmp,indsym,natom,0,nsym,qphon,dyfrlo_tmp2,symq,symrec) 226 dyfrxc(:,:,:)=dyfrlo_tmp2(1,:,:,:,1) 227 end if 228 ABI_DEALLOCATE(dyfrsym_tmp) 229 ABI_DEALLOCATE(dyfrlo_tmp2) 230 231 !Symmetrize nonlocal part of the dynamical matrix dyfrnl: 232 !atindx1 is used to reestablish the correct order of atoms 233 if (dyfr_nondiag==0) then 234 ABI_ALLOCATE(dyfrsym_tmp,(dyfr_cplex,3,3,natom,1)) 235 do iatom=1,natom 236 dyfrsym_tmp(:,:,:,atindx1(iatom),1)=dyfrnl(:,:,:,iatom,1) 237 end do 238 else 239 ABI_ALLOCATE(dyfrsym_tmp,(dyfr_cplex,3,3,natom,natom)) 240 do jatom=1,natom 241 do iatom=1,natom 242 dyfrsym_tmp(:,:,:,atindx1(iatom),atindx1(jatom))=dyfrnl(:,:,:,iatom,jatom) 243 end do 244 end do 245 end if 246 call dfpt_sydy(dyfr_cplex,dyfrsym_tmp,indsym,natom,dyfr_nondiag,nsym,qphon,dyfrnl,symq,symrec) 247 ABI_DEALLOCATE(dyfrsym_tmp) 248 249 !Collect local, nl xc core, and non-local part 250 !of the frozen wf dynamical matrix. 251 dyfrwf(:,:,:,:,:)=dyfrnl(:,:,:,:,:) 252 if (dyfr_nondiag==0) then 253 dyfrwf(1,:,:,:,1)=dyfrwf(1,:,:,:,1)+dyfrlo(:,:,:)+dyfrxc(:,:,:) 254 else 255 do iatom=1,natom 256 dyfrwf(1,:,:,iatom,iatom)=dyfrwf(1,:,:,iatom,iatom)+dyfrlo(:,:,iatom)+dyfrxc(:,:,iatom) 257 end do 258 end if 259 260 end subroutine dfpt_dyfro