TABLE OF CONTENTS


ABINIT/m_paw_pwaves_lmn [ Modules ]

[ Top ] [ Modules ]

NAME

  m_paw_pwaves_lmn

FUNCTION

  This module contains the definition of the paw_pwaves_lmn structured datatype,
  as well as related functions and methods.
  paw_pwaves_lmn is used to store the 3D values of the all-electron and of
  the pseudized part of the PAW partial waves on the set of FFT points falling
  inside the spheres around each atom.

COPYRIGHT

 Copyright (C) 2008-2018 ABINIT group (MG,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 .

NOTES

 * Routines tagged with "@type_name" are strongly connected to the definition of the data type.
   Strongly connected means that the proper functioning of the implementation relies on the
   assumption that the tagged procedure is consistent with the type declaration.
   Every time a developer changes the structure "type_name" adding new entries, he/she has to make sure
   that all the strongly connected routines are changed accordingly to accommodate the modification of the data type.
   Typical examples of strongly connected routines are creation, destruction or reset methods.

PARENTS

CHILDREN

SOURCE

33 #if defined HAVE_CONFIG_H
34 #include "config.h"
35 #endif
36 
37 #include "abi_common.h"
38 
39 MODULE m_paw_pwaves_lmn
40 
41  use defs_basis
42  use m_errors
43  use m_abicore
44  use m_xmpi
45  use m_sort
46 
47  use m_numeric_tools, only : wrap2_zero_one
48  use m_geometry,      only : xcart2xred
49  use m_pawrad,        only : pawrad_type, nderiv_gen, pawrad_deducer0
50  use m_pawtab,        only : pawtab_type
51  use m_pawfgrtab,     only : pawfgrtab_type
52  use m_paw_numeric,   only : paw_spline, paw_splint
53  use m_paw_sphharm,   only : initylmr
54  use m_paral_atom,    only : get_my_atmtab, free_my_atmtab
55 
56  implicit none
57 
58  private

m_paw_pwaves_lmn/paw_pwaves_lmn_free [ Functions ]

[ Top ] [ m_paw_pwaves_lmn ] [ Functions ]

NAME

 paw_pwaves_lmn_free

FUNCTION

INPUTS

OUTPUT

PARENTS

      classify_bands,exc_plot,m_paw_pwaves_lmn,m_wfd,pawmkaewf,screening
      sigma,wfk_analyze

CHILDREN

SOURCE

452 subroutine paw_pwaves_lmn_free(Paw_onsite)
453 
454 
455 !This section has been created automatically by the script Abilint (TD).
456 !Do not modify the following lines by hand.
457 #undef ABI_FUNC
458 #define ABI_FUNC 'paw_pwaves_lmn_free'
459 !End of the abilint section
460 
461  implicit none
462 
463 !Arguments ------------------------------------
464 !arrays
465  type(paw_pwaves_lmn_t),intent(inout) :: Paw_onsite(:)
466 
467 !Local variables-------------------------------
468 !scalars
469  integer :: iatom
470 
471 ! *************************************************************************
472 
473 !@paw_pwaves_lmn_t
474 
475  do iatom=LBOUND(Paw_onsite,DIM=1),UBOUND(Paw_onsite,DIM=1)
476    if (allocated(Paw_onsite(iatom)%phi)) then
477      ABI_FREE(Paw_onsite(iatom)%phi)
478    end if
479    if (allocated(Paw_onsite(iatom)%tphi)) then
480      ABI_FREE(Paw_onsite(iatom)%tphi)
481    end if
482    if (allocated(Paw_onsite(iatom)%r0shift)) then
483      ABI_FREE(Paw_onsite(iatom)%r0shift)
484    end if
485    if (allocated(Paw_onsite(iatom)%phi_gr )) then
486      ABI_FREE(Paw_onsite(iatom)%phi_gr)
487    end if
488    if (allocated(Paw_onsite(iatom)%tphi_gr)) then
489      ABI_FREE(Paw_onsite(iatom)%tphi_gr)
490    end if
491  end do
492 
493 end subroutine paw_pwaves_lmn_free

m_paw_pwaves_lmn/paw_pwaves_lmn_init [ Functions ]

[ Top ] [ m_paw_pwaves_lmn ] [ Functions ]

NAME

 paw_pwaves_lmn_init

FUNCTION

INPUTS

  mpi_atmtab(:)=--optional-- indexes of the atoms treated by current proc
  comm_atom=--optional-- MPI communicator over atoms

OUTPUT

PARENTS

      classify_bands,exc_plot,m_wfd,pawmkaewf,screening,sigma,wfk_analyze

CHILDREN

SOURCE

136 subroutine paw_pwaves_lmn_init(Paw_onsite,my_natom,natom,ntypat,rprimd,xcart,Pawtab, &
137 &                              Pawrad,local_pawfgrtab,optgrad,&
138 &                              mpi_atmtab,comm_atom) ! optional arguments (parallelism)
139 
140 
141 !This section has been created automatically by the script Abilint (TD).
142 !Do not modify the following lines by hand.
143 #undef ABI_FUNC
144 #define ABI_FUNC 'paw_pwaves_lmn_init'
145 !End of the abilint section
146 
147  implicit none
148 
149 !Arguments ------------------------------------
150 !scalars
151  integer,intent(in) :: my_natom,natom,ntypat
152  integer,optional,intent(in) :: optgrad,comm_atom
153 !arrays
154  integer,optional,target,intent(in) :: mpi_atmtab(:)
155  real(dp),intent(in) :: xcart(3,natom),rprimd(3,3)
156  type(pawtab_type),target,intent(in) :: Pawtab(ntypat)
157  type(pawrad_type),intent(in) :: Pawrad(ntypat)
158  type(pawfgrtab_type),intent(in) :: local_pawfgrtab(my_natom)
159  type(paw_pwaves_lmn_t),intent(out) :: Paw_onsite(my_natom)
160 
161 !Local variables-------------------------------
162 !scalars
163  integer :: itypat,ln_size,lmn_size,mesh_size,inl,iatom,iatom1,my_comm_atom,my_optgrad
164  integer :: nfgd,ifgd,ipsang,option_ylmr,normchoice,ii,jlmn,jl,jm,jlm,jln
165  logical :: my_atmtab_allocated,paral_atom
166  real(dp) :: phj,rR,tphj,ybcbeg,ybcend
167 !arrays
168  integer, allocatable :: iperm(:)
169  integer,ABI_CONTIGUOUS pointer :: indlmn(:,:)
170  integer,pointer :: my_atmtab(:)
171  real(dp) :: yvals(4),red(3),shift(3)
172  real(dp),allocatable :: ff(:),nrm(:),nrm_sort(:),phigrd(:,:),tphigrd(:,:),ylm_tmp(:,:),ylm(:,:),ylm_gr(:,:,:)
173  real(dp),allocatable :: rsph_red(:,:),rsph_cart(:,:),phigrd_gr(:,:),tphigrd_gr(:,:),gg(:)
174  type(paw_pwaves_lmn_t),allocatable :: Paw_lmn_spline(:)
175 
176 ! *************************************************************************
177 
178 !@paw_pwaves_lmn_t
179 
180  if (my_natom==0) return
181  my_optgrad = -1; if (present(optgrad)) my_optgrad = optgrad
182 
183  ABI_CHECK(all(local_pawfgrtab(:)%rfgd_allocated==1),"R vectors not allocated in pawfgrtab")
184 
185  ! Set up parallelism over atoms
186  paral_atom=(present(comm_atom).and.(my_natom/=natom))
187  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
188  my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
189  call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
190 
191  ! Prepare the spline. Calculate 2nd derivatives of partial waves for each atom type.
192  ABI_DT_MALLOC(Paw_lmn_spline,(ntypat))
193 
194  do itypat=1,ntypat
195    ln_size  =Pawtab(itypat)%basis_size
196    mesh_size=Pawtab(itypat)%mesh_size
197 
198    ABI_MALLOC(Paw_lmn_spline(itypat)%phi ,(mesh_size,ln_size))
199    ABI_MALLOC(Paw_lmn_spline(itypat)%tphi,(mesh_size,ln_size))
200 
201    do inl=1,ln_size ! Calculate 2nd derivatives of %phi and %tphi for each ln component.
202      ybcbeg=zero; ybcend=zero
203      call paw_spline(Pawrad(itypat)%rad,Pawtab(itypat)%phi(:,inl), mesh_size,&
204 &                    ybcbeg,ybcend,Paw_lmn_spline(itypat)%phi(:,inl))
205 
206      ybcbeg=zero; ybcend=zero
207      call paw_spline(Pawrad(itypat)%rad,Pawtab(itypat)%tphi(:,inl),mesh_size,&
208 &                    ybcbeg,ybcend,Paw_lmn_spline(itypat)%tphi(:,inl))
209    end do
210  end do
211  !
212  ! === spline for each atom ===
213  ! * FFT points within PAW sphere depend on the atom site.
214  do iatom1=1,my_natom
215    iatom=iatom1;if (paral_atom) iatom=my_atmtab(iatom1)
216 
217    itypat   = local_pawfgrtab(iatom1)%itypat
218    ln_size  = Pawtab(itypat)%basis_size
219    lmn_size = Pawtab(itypat)%lmn_size
220    mesh_size= Pawrad(itypat)%mesh_size
221    nfgd     = local_pawfgrtab(iatom1)%nfgd ! no. of points in the fine grid for this PAW sphere
222    indlmn   => Pawtab(itypat)%indlmn
223 
224    ! The points in the PAW sphere might belong to a different unit cell, in this case one has to
225    ! reconstruct the contribution to the AE psi_k given by the lattice-symmetric atom in another sphere.
226    ! Here I wrap rr back into the first unit cell keeping trace of the lattice vector L
227    ! needed so that rr = rr_first_ucell + L
228    ! The contribution to the AE u(r) in the first unit cell has to be multiplied by e^{-ikL}.
229 
230    Paw_onsite(iatom1)%nfgd     = nfgd
231    Paw_onsite(iatom1)%lmn_size = lmn_size
232    !
233    ABI_MALLOC(Paw_onsite(iatom1)%r0shift,(3,nfgd))
234 
235    ABI_MALLOC(rsph_red,(3,nfgd))
236    ABI_MALLOC(rsph_cart,(3,nfgd))
237    do ifgd=1,nfgd
238      rsph_cart(:,ifgd) = local_pawfgrtab(iatom1)%rfgd(:,ifgd) + xcart(:,iatom)
239    end do
240    call xcart2xred(nfgd,rprimd,rsph_cart,rsph_red) ! go to reduced coordinates.
241    ABI_FREE(rsph_cart)
242 
243    do ifgd=1,nfgd
244      call wrap2_zero_one(rsph_red(1,ifgd),red(1),shift(1)) ! rr = r_cell + shift
245      call wrap2_zero_one(rsph_red(2,ifgd),red(2),shift(2))
246      call wrap2_zero_one(rsph_red(3,ifgd),red(3),shift(3))
247      Paw_onsite(iatom1)%r0shift(:,ifgd) = NINT(shift)
248      !if (ANY( ABS(shift) > tol12)) then
249        !MSG_WARNING("rmR_red is outside the first unit cell.")
250        !write(ab_out,*)rsph_red(:,ifgd),shift
251      !end if
252    end do
253    ABI_FREE(rsph_red)
254    !
255    ! * Obtain |r-R| on fine grid, note that rfgd is given in Cartesian coordinates.
256    ABI_MALLOC(nrm,(nfgd))
257    do ifgd=1,nfgd
258      nrm(ifgd) = sqrt(dot_product(local_pawfgrtab(iatom1)%rfgd(:,ifgd),local_pawfgrtab(iatom1)%rfgd(:,ifgd)))
259    end do
260    !
261    ! * Compute Ylm for each r-R vector.
262    ipsang = 1 + (Pawtab(itypat)%l_size-1)/2 ! recall l_size=2*l_max-1 where l_max is shifted by 1.
263    ABI_MALLOC(ylm_tmp,(ipsang**2,nfgd))
264    normchoice = 1 ! Use computed norms of input vectors.
265    if (my_optgrad==1) then
266      option_ylmr=2 ! Compute Ylm(r-R) and its gradient
267      ABI_MALLOC(ylm_gr,(3,ipsang**2,nfgd))
268    else
269      option_ylmr= 1 ! To compute Ylm(r-R).
270      ABI_MALLOC(ylm_gr,(3,3,0))
271    end if
272    call initylmr(ipsang,normchoice,nfgd,nrm,option_ylmr,local_pawfgrtab(iatom1)%rfgd,ylm_tmp,ylm_gr)
273    !
274    !  Exchange dimensions for better memory access.
275    ABI_MALLOC(ylm,(nfgd,ipsang**2))
276    do ii=1,ipsang**2
277      ylm(:,ii) = ylm_tmp(ii,:)
278    end do
279    ABI_FREE(ylm_tmp)
280    !
281    ! In order to do spline fits, the |r-R| data must be sorted
282    ! Here we sort the nrm points, keeping track of which goes where
283    ABI_MALLOC(nrm_sort,(nfgd))
284    nrm_sort = nrm
285 
286    ABI_MALLOC(iperm,(nfgd))
287    do ifgd=1,nfgd
288      iperm(ifgd)=ifgd
289    end do
290 
291    call sort_dp(nfgd,nrm_sort,iperm,tol8)
292 
293    ! Now make spline fits of phi and tphi  onto the fine grid around the atom
294    ABI_MALLOC(phigrd,(nfgd,ln_size))
295    ABI_MALLOC(tphigrd,(nfgd,ln_size))
296    ABI_MALLOC(ff,(nfgd))
297    if (my_optgrad==1) then
298      ABI_MALLOC(phigrd_gr,(nfgd,ln_size))
299      ABI_MALLOC(tphigrd_gr,(nfgd,ln_size))
300      ABI_MALLOC(gg,(mesh_size))
301    end if
302 
303    do inl=1,ln_size
304      !
305      ! * splint phi onto points and reorder indices.
306      call paw_splint(mesh_size,Pawrad(itypat)%rad,Pawtab(itypat)%phi(:,inl),&
307 &                    Paw_lmn_spline(itypat)%phi(:,inl),nfgd,nrm_sort,ff)
308      do ifgd=1,nfgd
309        ii=iperm(ifgd)
310        phigrd(ii,inl) = ff(ifgd)
311      end do
312      !
313      ! * compute d phi/dr, interpolate onto points and reorder indices.
314      if (my_optgrad==1) then
315        ybcbeg=zero; ybcend=zero
316        call nderiv_gen(gg,Pawtab(itypat)%phi(:,inl),Pawrad(itypat))
317        call paw_spline(Pawrad(itypat)%rad,gg,mesh_size,ybcbeg,ybcend,&
318 &                      Paw_lmn_spline(itypat)%phi(:,inl))
319        call paw_splint(mesh_size,Pawrad(itypat)%rad,Paw_lmn_spline(itypat)%phi(:,inl),&
320 &                      Paw_lmn_spline(itypat)%phi(:,inl),nfgd,nrm_sort,ff)
321        do ifgd=1,nfgd
322          ii=iperm(ifgd)
323          phigrd_gr(ii,inl)  = ff(ifgd)
324        end do
325      end if
326      !
327      ! * compute d tphi/dr, interpolate onto points and reorder indices.
328      call paw_splint(mesh_size,Pawrad(itypat)%rad,Pawtab(itypat)%tphi(:,inl),&
329 &                    Paw_lmn_spline(itypat)%tphi(:,inl),nfgd,nrm_sort,ff)
330      do ifgd=1,nfgd
331        ii=iperm(ifgd)
332        tphigrd(ii,inl) = ff(ifgd)
333      end do
334      if (my_optgrad==1) then
335        ybcbeg=zero; ybcend=zero
336        call nderiv_gen(gg,Pawtab(itypat)%tphi(:,inl),Pawrad(itypat))
337        call paw_spline(Pawrad(itypat)%rad,gg,mesh_size,ybcbeg,ybcend,&
338 &                      Paw_lmn_spline(itypat)%tphi(:,inl))
339        call paw_splint(mesh_size,Pawrad(itypat)%rad,Paw_lmn_spline(itypat)%tphi(:,inl),&
340 &                      Paw_lmn_spline(itypat)%phi(:,inl),nfgd,nrm_sort,ff)
341        do ifgd=1,nfgd
342          ii=iperm(ifgd)
343          tphigrd_gr(ii,inl)  = ff(ifgd)
344        end do
345      end if
346    end do !inl
347 
348    ABI_FREE(ff)
349    if (my_optgrad==1) then
350      ABI_FREE(gg)
351    end if
352    !
353    ! === Calculate AE and PS partial waves inside the sphere ===
354    ! * recall that <r|phi>=u(r)*Slm(r^)/r, hence avoid division by zero except for s-waves.
355    ABI_MALLOC(Paw_onsite(iatom1)%phi ,(nfgd,lmn_size))
356    ABI_MALLOC(Paw_onsite(iatom1)%tphi,(nfgd,lmn_size))
357 
358    if (my_optgrad==1) then
359      ABI_MALLOC(Paw_onsite(iatom1)%phi_gr ,(3,nfgd,lmn_size))
360      ABI_MALLOC(Paw_onsite(iatom1)%tphi_gr,(3,nfgd,lmn_size))
361    end if
362 
363    do jlmn=1,lmn_size
364      jl  = indlmn(1,jlmn)
365      jm  = indlmn(2,jlmn)
366      jlm = indlmn(4,jlmn)
367      jln = indlmn(5,jlmn)
368 
369      do ifgd=1,nfgd ! loop over fine grid points in current PAW sphere
370        !if (nrm(ifgd)>tol16) then
371        if (nrm(ifgd)>tol10) then ! tol10 to be consistent with initylmr.
372          rR  = nrm(ifgd) ! value of |r-R|
373          !write(ab_out,*) 'rR:',rR,' phigrd:',phigrd(ifgd,jln),' tphigrd:',tphigrd(ifgd,jln),' ylm:',ylm(ifgd,jlm)
374          phj = phigrd (ifgd,jln)*ylm(ifgd,jlm)/rR
375          tphj= tphigrd(ifgd,jln)*ylm(ifgd,jlm)/rR
376          Paw_onsite(iatom1)%phi (ifgd,jlmn) =  phj
377          Paw_onsite(iatom1)%tphi(ifgd,jlmn) = tphj
378 
379          if (my_optgrad==1) then
380            Paw_onsite(iatom1)%phi_gr (1:3,ifgd,jlmn) = phigrd (ifgd,jln)*ylm_gr(1:3,jlm,ifgd) &
381 &            + phigrd_gr(ifgd,jln)*local_pawfgrtab(iatom1)%rfgd(1:3,ifgd)*ylm(ifgd,jlm)
382            Paw_onsite(iatom1)%tphi_gr (1:3,ifgd,jlmn) = tphigrd (ifgd,jln)*ylm_gr(1:3,jlm,ifgd) &
383 &            + tphigrd_gr(ifgd,jln)*local_pawfgrtab(iatom1)%rfgd(1:3,ifgd)*ylm(ifgd,jlm)
384          end if
385 
386        else ! Extrapolate if the point is at the origin
387          yvals(1) = zero
388          if (jl==0) then
389            yvals(2:4) = Pawtab(itypat)%phi(2:4,jln)/Pawrad(itypat)%rad(2:4)
390            call pawrad_deducer0(yvals,4,Pawrad(itypat))
391          end if
392          Paw_onsite(iatom1)%phi(ifgd,jlmn) = yvals(1) * ylm(ifgd,jlm)
393          yvals(1) = zero
394          if (jl==0) then
395            yvals(2:4) = Pawtab(itypat)%tphi(2:4,jln)/Pawrad(itypat)%rad(2:4)
396            call pawrad_deducer0(yvals,4,pawrad(itypat))
397          end if
398          Paw_onsite(iatom1)%tphi(ifgd,jlmn) = yvals(1) * ylm(ifgd,jlm)
399          ! The gradient is expected to go to zero at the origin
400          if (my_optgrad==1) then
401            Paw_onsite(iatom1)%phi_gr (1:3,ifgd,jlmn) = zero
402            Paw_onsite(iatom1)%tphi_gr (1:3,ifgd,jlmn) = zero
403          end if
404        end if
405 
406      end do !nfgd
407    end do !jlmn
408 
409    ABI_FREE(nrm)
410    ABI_FREE(nrm_sort)
411    ABI_FREE(iperm)
412    ABI_FREE(phigrd)
413    ABI_FREE(tphigrd)
414    ABI_FREE(ylm)
415    ABI_FREE(ylm_gr)
416    if (my_optgrad==1) then
417      ABI_FREE(phigrd_gr)
418      ABI_FREE(tphigrd_gr)
419    end if
420  end do !iatom
421  !
422  !* Free 2nd derivates used for spline.
423  call paw_pwaves_lmn_free(Paw_lmn_spline)
424  ABI_DATATYPE_DEALLOCATE(Paw_lmn_spline)
425 
426  ! Destroy atom table used for parallelism
427  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
428 
429 end subroutine paw_pwaves_lmn_init

m_paw_pwaves_lmn/paw_pwaves_lmn_t [ Types ]

[ Top ] [ m_paw_pwaves_lmn ] [ Types ]

NAME

 paw_pwaves_lmn_t

FUNCTION

  Datatype used to store the 3D values of the all-electron and of the pseudized part of the
  PAW partial waves on the set of FFT points falling inside the spheres around each atom.
  The data is mainly used for plotting the true PAW wavefunctions in real space.

SOURCE

 74  type,public ::  paw_pwaves_lmn_t
 75 
 76   integer :: nfgd
 77 
 78   integer :: lmn_size
 79 
 80   !$integer :: ngfft(18)
 81 
 82   integer,allocatable :: r0shift(:,:)
 83   ! r0shift(3,nfgd)
 84 
 85   !real(dp),allocatable:: phk_atm(:,:)
 86   ! phk_atmt(2,nfgd)
 87 
 88   real(dp),allocatable :: phi(:,:)
 89   ! phi (nfgd,lmn_size)
 90   ! \phi_{nlm}(ifgd) for each point of the FFT mesh located inside the PAW sphere (see pawfgrtab_type).
 91 
 92   real(dp),allocatable :: tphi(:,:)
 93   ! tphi (nfgd,lmn_size)
 94   ! \tphi_{nlm}(ifgd) for each point of the FFT mesh located inside the PAW sphere (see pawfgrtab_type).
 95 
 96   real(dp),allocatable :: phi_gr(:,:,:)
 97   ! phi_gr (3,nfgd,lmn_size)
 98   ! gradient, in cartesian coordinates, of \phi_{nlm}(ifgd) for each point of the FFT mesh
 99   ! located inside the PAW sphere (see pawfgrtab_type).
100 
101   real(dp),allocatable :: tphi_gr(:,:,:)
102   ! tphi_gr (3,nfgd,lmn_size)
103   ! gradient, in cartesian coordinates, of \tphi_{nlm}(ifgd) for each point of the FFT mesh
104   ! located inside the PAW sphere (see pawfgrtab_type).
105 
106  end type paw_pwaves_lmn_t
107 
108  public :: paw_pwaves_lmn_init
109  public :: paw_pwaves_lmn_free