TABLE OF CONTENTS
ABINIT/psp6cc_drh [ Functions ]
NAME
psp6cc_drh
FUNCTION
Compute the core charge density, for use in the XC core correction, following the function definition valid for the format 6 of pseudopotentials. Version modified by DHamann, with consistent treatment of the derivatives in this routine and the remaining of the code.
COPYRIGHT
Copyright (C) 1998-2017 ABINIT group (AF,DRH) 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
mmax=maximum number of points in real space grid in the psp file n1xccc=dimension of xccc1d ; 0 if no XC core correction is used rchrg=cut-off radius for the core density
OUTPUT
xccc1d(n1xccc,6)= 1D core charge function and its five first derivatives
PARENTS
psp6in
CHILDREN
cc_derivatives
NOTES
Test version by DRH - requires very smooth model core charge
SOURCE
39 #if defined HAVE_CONFIG_H 40 #include "config.h" 41 #endif 42 43 #include "abi_common.h" 44 45 subroutine psp6cc_drh(mmax,n1xccc,rchrg,xccc1d) 46 47 use defs_basis 48 use m_profiling_abi 49 use m_errors 50 51 !This section has been created automatically by the script Abilint (TD). 52 !Do not modify the following lines by hand. 53 #undef ABI_FUNC 54 #define ABI_FUNC 'psp6cc_drh' 55 use interfaces_64_psp, except_this_one => psp6cc_drh 56 !End of the abilint section 57 58 implicit none 59 60 !Arguments ------------------------------------ 61 !scalars 62 integer,intent(in) :: mmax,n1xccc 63 real(dp),intent(in) :: rchrg 64 !arrays 65 real(dp),intent(inout) :: xccc1d(n1xccc,6) !vz_i 66 67 !Local variables------------------------------- 68 !scalars 69 integer :: irad 70 character(len=500) :: errmsg 71 !arrays 72 real(dp),allocatable :: ff(:),ff1(:),ff2(:),rad(:) 73 74 !********************************************************************** 75 76 ABI_ALLOCATE(ff,(mmax)) 77 ABI_ALLOCATE(ff1,(mmax)) 78 ABI_ALLOCATE(ff2,(mmax)) 79 ABI_ALLOCATE(rad,(mmax)) 80 81 ! 82 !read from pp file the model core charge (ff) and first (ff1) and 83 !second (ff2) derivative on logarithmic mesh mmax; rad is the radial grid 84 !the input functions contain the 4pi factor, it must be rescaled. 85 86 !***drh test 87 write(std_out,'(a,2i6)') 'drh:psp6cc_drh - mmax,n1xccc',mmax,n1xccc 88 !***end drh test 89 do irad=1,mmax 90 read(tmp_unit,*,err=10,iomsg=errmsg) rad(irad),ff(irad),ff1(irad),ff2(irad) 91 ff(irad)=ff(irad)/4.d0/pi 92 ff1(irad)=ff1(irad)/4.d0/pi 93 ff2(irad)=ff2(irad)/4.d0/pi 94 end do 95 rad(1)=0.d0 96 97 call cc_derivatives(rad,ff,ff1,ff2,mmax,n1xccc,rchrg,xccc1d) 98 99 ABI_DEALLOCATE(ff) 100 ABI_DEALLOCATE(ff1) 101 ABI_DEALLOCATE(ff2) 102 ABI_DEALLOCATE(rad) 103 104 return 105 106 ! Handle IO error 107 10 continue 108 MSG_ERROR(errmsg) 109 110 end subroutine psp6cc_drh