TABLE OF CONTENTS
ABINIT/fresidrsp [ Functions ]
NAME
fresidrsp
FUNCTION
Compute the forces due to the residual of the potential (or density) in RECIPROCAL SPACE, using - the atomic density read in psp file (PAW or NC with nctval_spl e.g. psp8 format) - a gaussian atomic density (norm-conserving psps if nctval_spl is not available)
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (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 dtset <type(dataset_type)>=all input variables in this dataset | densty(ntypat,4)=parameters for initialisation of the density of each atom type | icoulomb=0 periodic treatment of Hartree potential, 1 use of Poisson solver | ixc= choice of exchange-correlation scheme | natom=number of atoms in cell. | nspden=number of spin-density components | typat(natom)=integer type for each atom in cell gmet(3,3)=reciprocal space metric gprimd(3,3)=reciprocal space dimensional primitive translations gsqcut=cutoff value on G**2 for sphere inside fft box mgfft=maximum size of 1D FFTs mpi_enreg=information about MPI parallelization mqgrid=number of grid pts in q array for atomic density spline n^AT(q) nattyp(ntypat)=number of atoms of each type in cell 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 ntypat=number of types of atoms in cell. 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)=1-dim phase information for given atom coordinates. qgrid(mqgrid)=q grid for spline atomic valence density n^AT(q) from 0 to qmax ucvol=unit cell volume (bohr**3). usepaw= 0 for non paw calculation; =1 for paw calculation vresid(nfft,nspden)=potential residual - (non-collinear magn. : only V11 and V22 are used) zion(ntypat)=charge on each type of atom (real number) znucl(ntypat)=atomic number, for each type of atom
OUTPUT
gresid(3,natom)=forces due to the residual of the potential
PARENTS
forces
CHILDREN
atm2fft,fourdp,wrtout
SOURCE
60 #if defined HAVE_CONFIG_H 61 #include "config.h" 62 #endif 63 64 #include "abi_common.h" 65 66 subroutine fresidrsp(atindx1,dtset,gmet,gprimd,gresid,gsqcut,mgfft,mpi_enreg,mqgrid,nattyp,nfft,& 67 & ngfft,ntypat,psps,pawtab,ph1d,qgrid,ucvol,usepaw,vresid,zion,znucl) 68 69 use defs_basis 70 use defs_abitypes 71 use m_profiling_abi 72 73 use defs_datatypes, only : pseudopotential_type 74 use m_atomdata, only : atom_length 75 use m_pawtab, only : pawtab_type 76 77 !This section has been created automatically by the script Abilint (TD). 78 !Do not modify the following lines by hand. 79 #undef ABI_FUNC 80 #define ABI_FUNC 'fresidrsp' 81 use interfaces_14_hidewrite 82 use interfaces_53_ffts 83 use interfaces_64_psp 84 !End of the abilint section 85 86 implicit none 87 88 !Arguments ------------------------------------ 89 !scalars 90 integer,intent(in) :: mgfft,mqgrid,nfft,ntypat,usepaw 91 real(dp),intent(in) :: gsqcut,ucvol 92 type(pseudopotential_type),intent(in) :: psps 93 type(MPI_type),intent(in) :: mpi_enreg 94 type(dataset_type),intent(in) :: dtset 95 !arrays 96 integer,intent(in) :: atindx1(dtset%natom),nattyp(ntypat),ngfft(18) 97 real(dp),intent(in) :: gmet(3,3),gprimd(3,3),ph1d(2,3*(2*mgfft+1)*dtset%natom) 98 real(dp),intent(in) :: qgrid(mqgrid),vresid(nfft,dtset%nspden),zion(ntypat) 99 real(dp),intent(in) :: znucl(ntypat) 100 real(dp),intent(out) :: gresid(3,dtset%natom) 101 type(pawtab_type),intent(in) :: pawtab(ntypat*usepaw) 102 103 !Local variables------------------------------- 104 !scalars 105 integer :: itypat,optatm,optdyfr,opteltfr,optgr,optn,optn2,optstr,optv 106 logical :: usegauss 107 !arrays 108 integer :: dummy3(3) 109 real(dp) :: dummy2(2) 110 real(dp) :: dummy_in1(0),dummy_in2(0) 111 real(dp) :: dummy_out1(0),dummy_out2(0),dummy_out3(0),dummy_out4(0),dummy_out5(0),dummy_out6(0) 112 real(dp) :: strn_dummy6(6),strv_dummy6(6) 113 real(dp),allocatable :: gauss(:,:),vresg(:,:),work(:) 114 115 ! ************************************************************************* 116 117 !Inits 118 optatm=0;optdyfr=0;opteltfr=0;optgr=1;optstr=0;optv=0;optn=1 119 ABI_ALLOCATE(vresg,(2,nfft)) 120 121 !Transfer potential residual to reciprocal space 122 !Use only Vres=Vres11+Vres22=Vres_up+Vres_dn 123 ABI_ALLOCATE(work,(nfft)) 124 work(:)=vresid(:,1) 125 if (dtset%nspden>=2) work(:)=work(:)+vresid(:,2) 126 call fourdp(1,vresg,work,-1,mpi_enreg,nfft,ngfft,dtset%paral_kgb,0) 127 ABI_DEALLOCATE(work) 128 129 !Determine wether a gaussan atomic density has to be used or not 130 usegauss=.true. 131 if (usepaw==0) usegauss = any(.not.psps%nctab(1:ntypat)%has_tvale) 132 if (usepaw==1) usegauss=(minval(pawtab(1:ntypat)%has_tvale)==0) 133 if (usegauss) then 134 optn2=3 135 ABI_ALLOCATE(gauss,(2,ntypat)) 136 do itypat=1,ntypat 137 gauss(1,itypat)=zion(itypat) 138 gauss(2,itypat) = atom_length(dtset%densty(itypat,1),zion(itypat),znucl(itypat)) 139 end do 140 call wrtout(std_out," Computing residual forces using gaussian functions as atomic densities", "COLL") 141 else 142 optn2=2 143 ABI_ALLOCATE(gauss,(2,0)) 144 call wrtout(std_out," Computing residual forces using atomic densities taken from pseudos", "COLL") 145 end if 146 147 !Compute forces due to residual 148 call atm2fft(atindx1,dummy_out1,dummy_out2,dummy_out3,dummy_out4,& 149 & dummy_out5,gauss,gmet,gprimd,gresid,dummy_out6,gsqcut,mgfft,& 150 & mqgrid,dtset%natom,nattyp,nfft,ngfft,ntypat,optatm,optdyfr,opteltfr,optgr,optn,optn2,optstr,optv,& 151 & psps,pawtab,ph1d,qgrid,dummy3,dummy_in1,strn_dummy6,strv_dummy6,ucvol,usepaw,vresg,vresg,vresg,dummy2,dummy_in2,& 152 & comm_fft=mpi_enreg%comm_fft,me_g0=mpi_enreg%me_g0,& 153 & paral_kgb=mpi_enreg%paral_kgb,distribfft=mpi_enreg%distribfft) 154 155 !In case of nspden>=2, has to apply 1/2 factor 156 if (dtset%nspden>=2) gresid=gresid*half 157 158 ABI_DEALLOCATE(gauss) 159 ABI_DEALLOCATE(vresg) 160 161 end subroutine fresidrsp