TABLE OF CONTENTS
ABINIT/upf2abinit [ Functions ]
NAME
upf2abinit
FUNCTION
This routine wraps a call to a PWSCF module, which reads in a UPF (PWSCF / Espresso) format pseudopotential, then transfers data to abinit internal variables. "UPF PWSCF format" (pspcod=11)
COPYRIGHT
Copyright (C) 2009-2018 ABINIT group (MJV) 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
filpsp = name of file with UPF data psps = sturcture with global dimension data for pseudopotentials, header info ...
OUTPUT
pspxc = index of xc functional for this pseudo lmax_ = maximal angular momentum lloc = local component chosen for pseudopotential mmax = maximum number of points in real space radial grid znucl = charge of species nucleus zion = valence charge epsatm = integral of local potential - coulomb potential of zion xcccrc = radius for non linear core correction ekb(dimekb)= Kleinman Bylander energies, see pspatm.F90 epsatm=$ (4\pi)\int_0^\infty [r^2 (V(r)+\frac{Zv}{r}) dr]$(hartree) indlmn(6,i)= array giving l,m,n,lm,ln,s for i=ln (if useylm=0) or i=lmn (if useylm=1) vlspl(mqgrid_vl,2)=q^2 Vloc(q) and second derivatives from spline fit ffspl(mqgrid_ff,2,lnmax)=Kleinman-Bylander form factor f_l(q) and second derivative from spline fit for each angular momentum and each projector; if any, spin-orbit components begin at l=mpsang+1 nproj= number of projectors for each channel xccc1d(n1xccc*(1-usepaw),6)=1D core charge function and five derivatives, from psp file (used in NC only)
SIDE EFFECTS
NOTES
PARENTS
pspatm
CHILDREN
atomdata_from_symbol,cc_derivatives,nderiv,psp11lo,psp11nl,read_pseudo smooth,spline,upfxc2abi
SOURCE
58 #if defined HAVE_CONFIG_H 59 #include "config.h" 60 #endif 61 62 #include "abi_common.h" 63 64 subroutine upf2abinit (filpsp, znucl, zion, pspxc, lmax_, lloc, mmax, & 65 & psps, epsatm, xcccrc, indlmn, ekb, ffspl, nproj_l, vlspl, xccc1d) 66 67 use defs_basis 68 use defs_datatypes 69 use m_splines 70 use m_profiling_abi 71 use m_errors 72 use m_atomdata 73 use pseudo_pwscf ! pwscf module with all data explicit! 74 75 use m_io_tools, only : open_file 76 use m_numeric_tools, only : smooth, nderiv 77 78 !This section has been created automatically by the script Abilint (TD). 79 !Do not modify the following lines by hand. 80 #undef ABI_FUNC 81 #define ABI_FUNC 'upf2abinit' 82 use interfaces_11_qespresso_ext 83 use interfaces_57_iopsp_parser 84 use interfaces_64_psp, except_this_one => upf2abinit 85 !End of the abilint section 86 87 implicit none 88 89 !Arguments ------------------------------- 90 91 character(len=fnlen), intent(in) :: filpsp 92 type(pseudopotential_type),intent(in) :: psps 93 ! used contents: 94 ! psps%lmnmax 95 ! psps%mqgrid_ff 96 ! psps%mqgrid_vl 97 ! psps%dimekb 98 ! psps%n1xccc 99 ! psps%qgrid_ff 100 ! psps%qgrid_vl 101 102 integer, intent(out) :: pspxc, lmax_, lloc, mmax 103 real(dp), intent(out) :: znucl, zion 104 real(dp), intent(out) :: epsatm, xcccrc 105 !arrays 106 integer, intent(out) :: indlmn(6,psps%lmnmax) 107 integer, intent(out) :: nproj_l(psps%mpssoang) 108 real(dp), intent(inout) :: ekb(psps%dimekb) !vz_i 109 real(dp), intent(inout) :: ffspl(psps%mqgrid_ff,2,psps%lnmax) !vz_i 110 real(dp), intent(out) :: vlspl(psps%mqgrid_vl,2) 111 real(dp), intent(inout) :: xccc1d(psps%n1xccc,6) !vz_i 112 113 !Local variables ------------------------- 114 integer :: ir, iproj, ll, iunit 115 real(dp) :: yp1, ypn 116 character(len=500) :: msg 117 type(atomdata_t) :: atom 118 119 logical, allocatable :: found_l(:) 120 real(dp), allocatable :: work_space(:),work_spl(:) 121 real(dp), allocatable :: ff(:), ff1(:), ff2(:), rad_cc(:), proj(:,:) 122 123 ! ######### in module pseudo: ############ 124 ! 125 ! only npsx = 1 is used here 126 ! grids are allocated for much larger fixed length (ndm=2000) 127 ! number of species (6) and projectors (8) as well... 128 ! 129 ! psd(npsx) = specied string 130 ! pseudotype = uspp / nc string 131 ! dft(npsx) = exchange correlation string (20 chars) 132 ! lmax(npsx) = maximum l channel 133 ! mesh(npsx) = number of points for local pot 134 ! nbeta(npsx) = number of projectors (beta functions for uspp) 135 ! nlcc(npsx) = flag for presence of NL core correction 136 ! zp(npsx) = valence ionic charge 137 ! r(ndm,npsx) = radial mesh 138 ! rab(ndm,npsx) = dr / di for radial mesh 139 ! rho_atc(ndm,npsx) = NLCC pseudocharge density 140 ! vloc0(ndm,npsx) = local pseudopotential 141 ! betar(ndm, nbrx, npsx) = projector functions in real space mesh 142 ! lll(nbrx,npsx) = angular momentum channel for each projector 143 ! ikk2(nbrx,npsx) = maximum index for each projector function 144 ! dion(nbrx,nbrx,npsx) = dij or Kleinman Bylander energies 145 ! 146 ! ######## end description of pseudo module contents ########## 147 148 ! ********************************************************************* 149 150 !call pwscf routine for reading in UPF 151 if (open_file (filpsp,msg,newunit=iunit,status='old',form='formatted') /= 0) then 152 MSG_ERROR(msg) 153 end if 154 155 !read in psp data to static data in pseudo module, for ipsx == 1 156 call read_pseudo(1,iunit) 157 close (iunit) 158 159 !convert to Ha units 160 vloc0 = half * vloc0 161 !betar = half * betar ! ??? 162 dion = half * dion 163 164 !if upf file is a USPP one, stop 165 if (pseudotype == 'US') then 166 MSG_ERROR('upf2abinit: USPP UPF files not supported') 167 end if 168 169 !copy over to abinit internal arrays and vars 170 call upfxc2abi(dft(1), pspxc) 171 lmax_ = lmax(1) 172 173 !Check if the local component is one of the angular momentum channels 174 !effectively if one of the ll is absent from the NL projectors 175 ABI_ALLOCATE(found_l,(0:lmax_)) 176 found_l = .true. 177 do ll = 0, lmax_ 178 if (any(lll(1:nbeta(1),1) == ll)) then 179 found_l(ll) = .false. 180 end if 181 end do 182 if (count(found_l) /= 1) then 183 lloc = -1 184 else 185 do ll = 0, lmax_ 186 if (found_l(ll)) then 187 lloc = ll 188 exit 189 end if 190 end do 191 end if 192 ABI_DEALLOCATE(found_l) 193 !FIXME: do something about lloc == -1 194 195 call atomdata_from_symbol(atom,psd(1)) 196 znucl = atom%znucl 197 zion = zp(1) 198 mmax = mesh(1) 199 200 call psp11lo(rab(1:mmax,1),epsatm,mmax,psps%mqgrid_vl,psps%qgrid_vl,& 201 & vlspl(:,1),r(1:mmax,1),vloc0(1:mmax,1),yp1,ypn,zion) 202 203 !Fit spline to q^2 V(q) (Numerical Recipes subroutine) 204 ABI_ALLOCATE(work_space,(psps%mqgrid_vl)) 205 ABI_ALLOCATE(work_spl,(psps%mqgrid_vl)) 206 call spline (psps%qgrid_vl,vlspl(:,1),psps%mqgrid_vl,yp1,ypn,work_spl) 207 vlspl(:,2)=work_spl(:) 208 ABI_DEALLOCATE(work_space) 209 ABI_DEALLOCATE(work_spl) 210 211 !this has to do the FT of the projectors to reciprocal space 212 ! allocate proj to avoid temporary copy. 213 ABI_ALLOCATE(proj, (mmax,1:nbeta(1))) 214 proj = betar(1:mmax,1:nbeta(1),1) 215 216 call psp11nl(ffspl, indlmn, mmax, psps%lnmax, psps%lmnmax, psps%mqgrid_ff, & 217 & nbeta(1), proj, lll(1:nbeta(1),1), ikk2(1:nbeta(1),1), & 218 & psps%qgrid_ff, r(1:mmax,1), rab(1:mmax,1), psps%useylm) 219 220 ABI_FREE(proj) 221 222 nproj_l = 0 223 do iproj = 1, nbeta(1) 224 ll = lll(iproj,1) 225 nproj_l(ll+1) = nproj_l(ll+1) + 1 226 end do 227 228 !shape = dimekb vs. shape = n_proj 229 do ll = 1, nbeta(1) 230 ekb(ll) = dion(ll,ll,1) 231 end do 232 233 xcccrc = zero 234 xccc1d = zero 235 !if we find a core density, do something about it 236 !rho_atc contains the nlcc density 237 !rho_at contains the total density 238 if (nlcc(1)) then 239 ABI_ALLOCATE(ff,(mmax)) 240 ABI_ALLOCATE(ff1,(mmax)) 241 ABI_ALLOCATE(ff2,(mmax)) 242 ff(1:mmax) = rho_atc(1:mmax,1) ! model core charge without derivative factor 243 244 ff1 = zero 245 call nderiv(one,ff,ff1,mmax,1) ! first derivative 246 ff1(1:mmax) = ff1(1:mmax) / rab(1:mmax,1) 247 call smooth(ff1, mmax, 15) ! run 15 iterations of smoothing 248 249 ff2 = zero 250 call nderiv(one, ff1, ff2, mmax, 1) ! second derivative 251 ff2(1:mmax) = ff2(1:mmax) / rab(1:mmax,1) 252 call smooth(ff2, mmax, 15) ! run 10 iterations of smoothing? 253 254 ! determine a good rchrg = xcccrc 255 do ir = mmax, 1, -1 256 if (abs(ff(ir)) > 1.e-6) then 257 xcccrc = r(ir,1) 258 exit 259 end if 260 end do 261 ABI_ALLOCATE(rad_cc,(mmax)) 262 rad_cc = r(1:mmax,1) 263 rad_cc(1) = zero ! force this so that the core charge covers whole spline interval. 264 call cc_derivatives(rad_cc,ff,ff1,ff2,mmax,psps%n1xccc,xcccrc,xccc1d) 265 ABI_DEALLOCATE(rad_cc) 266 ABI_DEALLOCATE(ff) 267 ABI_DEALLOCATE(ff1) 268 ABI_DEALLOCATE(ff2) 269 270 end if !if nlcc present 271 272 end subroutine upf2abinit