TABLE OF CONTENTS
ABINIT/pawgylmg [ Functions ]
NAME
pawgylmg
FUNCTION
PAW: Compute Fourier transform of each g_l(r).Y_lm(r) function
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 .
INPUTS
gprimd(3,3)=dimensional reciprocal space primitive translations kg(3,npw)=integer coordinates of planewaves in basis sphere for this k point. kpg(npw,nkpg)= (k+G) components (only if useylm=1) kpt(3)=reduced coordinates of k point lmax=1+max. value of l angular momentum nkpg=second dimension of kpg_k (0 if useylm=0) npw=number of planewaves in basis sphere ntypat=number of types of atoms pawtab(ntypat) <type(pawtab_type)>=paw tabulated starting data ylm(npw,lmax**2)=real spherical harmonics for each G and k point
OUTPUT
gylmg(npw,lmax**2,ntypat)=Fourier transform of each g_l(r).Y_lm(r) function
PARENTS
suscep_stat
CHILDREN
splfit
SOURCE
38 #if defined HAVE_CONFIG_H 39 #include "config.h" 40 #endif 41 42 #include "abi_common.h" 43 44 subroutine pawgylmg(gprimd,gylmg,kg,kpg,kpt,lmax,nkpg,npw,ntypat,pawtab,ylm) 45 46 use m_profiling_abi 47 48 use defs_basis 49 use m_errors 50 use m_splines 51 use m_pawtab, only : pawtab_type 52 53 !This section has been created automatically by the script Abilint (TD). 54 !Do not modify the following lines by hand. 55 #undef ABI_FUNC 56 #define ABI_FUNC 'pawgylmg' 57 !End of the abilint section 58 59 implicit none 60 61 !Arguments ------------------------------------ 62 !scalars 63 integer,intent(in) :: lmax,nkpg,npw,ntypat 64 !arrays 65 integer,intent(in) :: kg(3,npw) 66 real(dp),intent(in) :: gprimd(3,3),kpg(npw,nkpg),kpt(3) 67 real(dp),intent(in) :: ylm(npw,lmax**2) 68 type(pawtab_type),intent(in) :: pawtab(ntypat) 69 70 real(dp),intent(out) :: gylmg(npw,lmax**2,ntypat) 71 72 !Local variables------------------------------- 73 !scalars 74 integer :: ig,ilm,itypat,ll,l0,mm,mqgrid 75 real(dp) :: kpg1,kpg2,kpg3,kpgc1,kpgc2,kpgc3 76 !arrays 77 real(dp),allocatable :: glg(:),qgrid(:),kpgnorm(:),shpf(:,:),work(:) 78 79 ! ************************************************************************* 80 81 DBG_ENTER("COLL") 82 83 !Get |k+G|: 84 ABI_ALLOCATE(kpgnorm,(npw)) 85 if (nkpg<3) then 86 do ig=1,npw 87 kpg1=kpt(1)+dble(kg(1,ig));kpg2=kpt(2)+dble(kg(2,ig));kpg3=kpt(3)+dble(kg(3,ig)) 88 kpgc1=kpg1*gprimd(1,1)+kpg2*gprimd(1,2)+kpg3*gprimd(1,3) 89 kpgc2=kpg1*gprimd(2,1)+kpg2*gprimd(2,2)+kpg3*gprimd(2,3) 90 kpgc3=kpg1*gprimd(3,1)+kpg2*gprimd(3,2)+kpg3*gprimd(3,3) 91 kpgnorm(ig)=sqrt(kpgc1*kpgc1+kpgc2*kpgc2+kpgc3*kpgc3) 92 end do 93 else 94 do ig=1,npw 95 kpgc1=kpg(ig,1)*gprimd(1,1)+kpg(ig,2)*gprimd(1,2)+kpg(ig,3)*gprimd(1,3) 96 kpgc2=kpg(ig,1)*gprimd(2,1)+kpg(ig,2)*gprimd(2,2)+kpg(ig,3)*gprimd(2,3) 97 kpgc3=kpg(ig,1)*gprimd(3,1)+kpg(ig,2)*gprimd(3,2)+kpg(ig,3)*gprimd(3,3) 98 kpgnorm(ig)=sqrt(kpgc1*kpgc1+kpgc2*kpgc2+kpgc3*kpgc3) 99 end do 100 end if 101 102 ABI_ALLOCATE(glg,(npw)) 103 ABI_ALLOCATE(work,(npw)) 104 105 write(std_out,*) ' lmax, pawtab(:)%lcut_size ', lmax, pawtab(:)%lcut_size 106 107 !Loop over types of atoms 108 do itypat=1,ntypat 109 110 mqgrid=pawtab(itypat)%mqgrid_shp 111 ABI_ALLOCATE(qgrid,(mqgrid)) 112 ABI_ALLOCATE(shpf,(mqgrid,2)) 113 qgrid(1:mqgrid)=pawtab(itypat)%qgrid_shp(1:mqgrid) 114 115 ! Loops over (l,m) values 116 do ll=0,pawtab(itypat)%lcut_size-1 117 l0=ll**2+ll+1 118 119 shpf(1:mqgrid,1:2)=pawtab(itypat)%shapefncg(1:mqgrid,1:2,1+ll) 120 call splfit(qgrid,work,shpf,0,kpgnorm,glg,mqgrid,npw) 121 122 do mm=-ll,ll 123 ilm=l0+mm 124 125 gylmg(1:npw,ilm,itypat)=ylm(1:npw,ilm)*glg(1:npw) 126 127 ! End loops over (l,m) values 128 end do 129 end do 130 131 ! End loop over atom types 132 ABI_DEALLOCATE(qgrid) 133 ABI_DEALLOCATE(shpf) 134 end do 135 136 ABI_DEALLOCATE(kpgnorm) 137 ABI_DEALLOCATE(glg) 138 ABI_DEALLOCATE(work) 139 140 DBG_EXIT("COLL") 141 142 end subroutine pawgylmg