TABLE OF CONTENTS


ABINIT/dfpt_dyfro [ Functions ]

[ Top ] [ 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