./65_paw/nhatgrid.F90

ABINIT/nhatgrid

ABINIT/nhatgrid [ Functions ]

NAME

 nhatgrid

FUNCTION

 Determine parts of the rectangular (fine) grid that are contained
 inside spheres around atoms (used to compute n_hat density).
 If corresponding option is selected, compute also g_l(r)*Y_lm(r)
 (and derivatives) on this grid (g_l=radial shape function).

 Copyright (C) 1998-2018 ABINIT group (FJ, 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 (see gstate.f)
  distribfft<type(distribfft_type)>=--optional-- contains all the informations related
                                    to the FFT parallelism and plane sharing
  gmet(3,3)=reciprocal space metric tensor in bohr**-2
  mpi_atmtab(:)=--optional-- indexes of the atoms treated by current proc
  comm_atom=--optional-- MPI communicator over atoms
  comm_fft=--optional-- MPI communicator over FFT components
  my_natom=number of atoms treated by current processor
  natom=total number of atoms in cell
  nattyp(ntypat)= # atoms of each type.
  ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft
  ntypat=number of types of atoms in unit cell
  optcut= option for the cut-off radius of spheres:
          if optcut=0, cut-off radius=pawtab%rshp=cut-off radius of compensation charge
          if optcut=1, cut-off radius=pawtab%rpaw=radius of PAW augmentation regions
  optgr0= 1 if g_l(r)*Y_lm(r) are computed
  optgr1= 1 if first derivatives of g_l(r)*Y_lm(r) are computed
  optgr2= 1 if second derivatives of g_l(r)*Y_lm(r) are computed
  optrad= 1 if vectors (r-r_atom) on the fine grid around atoms have to be stored
  pawfgrtab(natom) <type(pawfgrtab_type)>=atomic data given on fine rectangular grid
  pawtab(ntypat) <type(pawtab_type)>=paw tabulated starting data
  rprimd(3,3)=dimensional primitive translations in real space (bohr)
  typat(natom)=type (integer) for each atom
  typord=1 if the output is ordered by type of atoms, 0 otherwise
  ucvol=unit cell volume in bohr**3
  xred(3,natom)=reduced dimensionless atomic coordinates

OUTPUT

  pawfgrtab(natom)%ifftsph(nfgd)=FFT index (fine grid) of a points in paw spheres around each atom
  pawfgrtab(natom)%nfgd= number of (fine grid) FFT points in paw spheres around atoms
  if (optgr0==1)
    pawfgrtab(natom)%gylm(nfgd,l_size**2)= g_l(r)*Y_lm(r) around each atom
  if (optgr1==1)
    pawfgrtab(natom)%gylmgr(3,nfgd,l_size**2)= derivatives of g_l(r)*Y_lm(r) wrt cart. coordinates
  if (optgr2==1)
    pawfgrtab(natom)%gylmgr2(6,nfgd,l_size**2)= second derivatives of g_l(r)*Y_lm(r) wrt cart. coordinates
  if (optrad==1)
    pawfgrtab(natom)%rfgd(3,nfgd)= coordinates of r-r_atom around each atom

PARENTS

      afterscfloop,bethe_salpeter,classify_bands,denfgr,exc_plot,m_wfd
      pawmkaewf,respfn,scfcv,screening,sigma,wfk_analyze

CHILDREN

      free_my_atmtab,get_my_atmtab,pawgylm,pawrfgd_fft,timab

SOURCE

 68 #if defined HAVE_CONFIG_H
 69 #include "config.h"
 70 #endif
 71 
 72 #include "abi_common.h"
 73 
 74 subroutine nhatgrid(atindx1,gmet,my_natom,natom,nattyp,ngfft,ntypat,&
 75 & optcut,optgr0,optgr1,optgr2,optrad,pawfgrtab,pawtab,rprimd,typat,ucvol,xred, &
 76 & mpi_atmtab,comm_atom,comm_fft,distribfft,typord) ! optional arguments (parallelism)
 77 
 78  use defs_basis
 79  use m_profiling_abi
 80  use m_errors
 81 
 82  use m_xmpi,         only : xmpi_comm_self,xmpi_comm_rank,xmpi_comm_size
 83  use m_pawtab,       only : pawtab_type
 84  use m_pawfgrtab,    only : pawfgrtab_type
 85  use m_paw_finegrid, only : pawgylm, pawrfgd_fft
 86  use m_paral_atom,   only : get_my_atmtab, free_my_atmtab
 87  use m_distribfft,   only : distribfft_type
 88 
 89 !This section has been created automatically by the script Abilint (TD).
 90 !Do not modify the following lines by hand.
 91 #undef ABI_FUNC
 92 #define ABI_FUNC 'nhatgrid'
 93  use interfaces_18_timing
 94 !End of the abilint section
 95 
 96  implicit none
 97 
 98 !Arguments ---------------------------------------------
 99 !scalars
100  integer,intent(in) :: my_natom,natom,ntypat,optcut,optgr0,optgr1,optgr2,optrad
101  integer,optional,intent(in) :: comm_atom,comm_fft,typord
102  real(dp),intent(in) :: ucvol
103  type(distribfft_type),optional,target,intent(in)  :: distribfft
104 !arrays
105  integer,intent(in) :: ngfft(18),typat(natom)
106  integer,intent(in),target :: atindx1(natom),nattyp(ntypat)
107  integer,optional,target,intent(in) :: mpi_atmtab(:)
108  real(dp),intent(in) :: gmet(3,3),rprimd(3,3),xred(3,natom)
109  type(pawfgrtab_type),intent(inout) :: pawfgrtab(my_natom)
110  type(pawtab_type),intent(in) :: pawtab(ntypat)
111 
112 !Local variables ------------------------------
113 !scalars
114  integer :: i3,iat,iatm,iatom,iatom_,iatom_tot,itypat,lm_size,me_fft,my_comm_atom,n1,n2,n3,nfgd
115  logical :: grid_found,my_atmtab_allocated,paral_atom
116  real(dp) :: rcut
117  character(len=500) :: msg
118 !arrays
119  integer,allocatable :: ifftsph_tmp(:)
120  integer,pointer :: my_atindx1(:),my_atmtab(:),my_nattyp(:)
121  integer, ABI_CONTIGUOUS pointer :: fftn3_distrib(:),ffti3_local(:)
122  real(dp) :: tsec(2)
123  real(dp),allocatable :: rfgd_tmp(:,:)
124 
125 ! *************************************************************************
126 
127  DBG_ENTER("COLL")
128 
129  call timab(559,1,tsec)
130  if (my_natom==0) return
131 
132 !Set up parallelism over FFT
133  me_fft=0
134  if (present(comm_fft)) then
135    me_fft=xmpi_comm_rank(comm_fft)
136  end if
137 
138 !Set up parallelism over atoms
139  paral_atom=(present(comm_atom).and.(my_natom/=natom))
140  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
141  my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
142  call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,&
143 & my_natom_ref=my_natom)
144  if (paral_atom) then
145    ABI_ALLOCATE(my_atindx1,(natom))
146    ABI_ALLOCATE(my_nattyp,(ntypat))
147    my_atindx1(:)=0;my_nattyp(:)=0
148    iat=1
149    do itypat=1,ntypat
150      if (my_natom>0) then
151        do iatom=1,my_natom
152          if(typat(my_atmtab(iatom))==itypat)then
153            my_nattyp(itypat)=my_nattyp(itypat)+1
154            my_atindx1(iat)=iatom
155            iat=iat+1
156          end if
157        end do
158      end if
159    end do
160  else
161    my_atindx1 => atindx1
162    my_nattyp => nattyp
163  end if
164 
165 !Get the distrib associated with this fft_grid
166  n1=ngfft(1);n2=ngfft(2);n3=ngfft(3)
167  if (present(distribfft)) then
168    grid_found=.false.
169    if (n2 == distribfft%n2_coarse) then
170      if (n3== size(distribfft%tab_fftdp3_distrib)) then
171        fftn3_distrib => distribfft%tab_fftdp3_distrib
172        ffti3_local => distribfft%tab_fftdp3_local
173        grid_found=.true.
174      end if
175    end if
176    if (n2 == distribfft%n2_fine) then
177      if (n3 == size(distribfft%tab_fftdp3dg_distrib)) then
178        fftn3_distrib => distribfft%tab_fftdp3dg_distrib
179        ffti3_local => distribfft%tab_fftdp3dg_local
180        grid_found = .true.
181      end if
182    end if
183    if (.not.(grid_found)) then
184      msg='Unable to find an allocated distrib for this fft grid!'
185      MSG_BUG(msg)
186    end if
187  else
188    ABI_ALLOCATE(fftn3_distrib,(n3))
189    ABI_ALLOCATE(ffti3_local,(n3))
190    fftn3_distrib=0;ffti3_local=(/(i3,i3=1,n3)/)
191  end if
192 
193 !Loop over types of atom
194 !-------------------------------------------
195  iatm=0
196  do itypat=1,ntypat
197 
198    if (optcut==1) then
199      rcut=pawtab(itypat)%rpaw
200    else
201      rcut=pawtab(itypat)%rshp
202    end if
203 
204 !  Loop over atoms
205 !  -------------------------------------------
206    do iat=1,my_nattyp(itypat)
207      iatm=iatm+1;iatom=my_atindx1(iatm)
208      iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
209      iatom_=iatom;if(present(typord)) iatom_=merge(iatm,iatom,typord==1)
210      lm_size=pawfgrtab(iatom_)%l_size**2
211 
212 !    ------------------------------------------------------------------
213 !    A-Determine FFT points and r-R vectors around the atom
214 !    ------------------------------------------------------------------
215 
216      call pawrfgd_fft(ifftsph_tmp,gmet,n1,n2,n3,nfgd,rcut,rfgd_tmp,rprimd,ucvol,&
217 &     xred(:,iatom_tot),fft_distrib=fftn3_distrib,fft_index=ffti3_local,me_fft=me_fft)
218 
219 !    Allocate arrays defining sphere (and related data) around current atom
220      if (allocated(pawfgrtab(iatom_)%ifftsph)) then
221        ABI_DEALLOCATE(pawfgrtab(iatom_)%ifftsph)
222      end if
223      ABI_ALLOCATE(pawfgrtab(iatom_)%ifftsph,(nfgd))
224      pawfgrtab(iatom_)%nfgd=nfgd
225      pawfgrtab(iatom_)%ifftsph(1:nfgd)=ifftsph_tmp(1:nfgd)
226 
227      if (optrad==1) then
228        if (allocated(pawfgrtab(iatom_)%rfgd))  then
229          ABI_DEALLOCATE(pawfgrtab(iatom_)%rfgd)
230        end if
231        ABI_ALLOCATE(pawfgrtab(iatom_)%rfgd,(3,nfgd))
232        pawfgrtab(iatom_)%rfgd_allocated=1
233        pawfgrtab(iatom_)%rfgd(1:3,1:nfgd)=rfgd_tmp(1:3,1:nfgd)
234      end if
235 
236      if (optgr0==1) then
237        if (allocated(pawfgrtab(iatom_)%gylm))  then
238          ABI_DEALLOCATE(pawfgrtab(iatom_)%gylm)
239        end if
240        ABI_ALLOCATE(pawfgrtab(iatom_)%gylm,(nfgd,lm_size))
241        pawfgrtab(iatom_)%gylm_allocated=1
242      end if
243 
244      if (optgr1==1) then
245        if (allocated(pawfgrtab(iatom_)%gylmgr))  then
246          ABI_DEALLOCATE(pawfgrtab(iatom_)%gylmgr)
247        end if
248        ABI_ALLOCATE(pawfgrtab(iatom_)%gylmgr,(3,nfgd,lm_size))
249        pawfgrtab(iatom_)%gylmgr_allocated=1
250      end if
251 
252      if (optgr2==1) then
253        if (allocated(pawfgrtab(iatom_)%gylmgr2))  then
254          ABI_DEALLOCATE(pawfgrtab(iatom_)%gylmgr2)
255        end if
256        ABI_ALLOCATE(pawfgrtab(iatom_)%gylmgr2,(6,nfgd,lm_size))
257        pawfgrtab(iatom_)%gylmgr2_allocated=1
258      end if
259 
260 !    ------------------------------------------------------------------
261 !    B-Calculate g_l(r-R)*Y_lm(r-R) for each r around the atom R
262 !    ------------------------------------------------------------------
263      if (optgr0+optgr1+optgr2>0) then
264        call pawgylm(pawfgrtab(iatom_)%gylm,pawfgrtab(iatom_)%gylmgr,pawfgrtab(iatom_)%gylmgr2,&
265 &       lm_size,nfgd,optgr0,optgr1,optgr2,pawtab(itypat),rfgd_tmp(:,1:nfgd))
266      end if
267 
268 !    End loops over types/atoms
269 !    -------------------------------------------
270      ABI_DEALLOCATE(ifftsph_tmp)
271      ABI_DEALLOCATE(rfgd_tmp)
272    end do
273  end do
274 
275 !Destroy atom tables used for parallelism
276  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
277  if (paral_atom) then
278    ABI_DEALLOCATE(my_atindx1)
279    ABI_DEALLOCATE(my_nattyp)
280  end if
281 
282  if (.not.present(distribfft)) then
283    ABI_DEALLOCATE(fftn3_distrib)
284    ABI_DEALLOCATE(ffti3_local)
285  end if
286 
287  call timab(559,2,tsec)
288 
289  DBG_EXIT("COLL")
290 
291 end subroutine nhatgrid

TABLE OF CONTENTS

ABINIT/nhatgrid [ Functions ]