TABLE OF CONTENTS
ABINIT/read_atomden [ Functions ]
NAME
read_atomden
FUNCTION
COPYRIGHT
Copyright (C) 2005-2018 ABINIT group (SM,VR,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 .
INPUTS
natom : number of atoms in cell nspden : number of spin densities ntypat : number of types of atoms in the cell pawfgr <type(pawfgr_type)> : data about the fine grid typat(natom) : list of atom types
OUTPUT
rhor_paw(pawfgr%nfft,nspden) : full electron density on the fine grid
SIDE EFFECTS
NOTES
PARENTS
outscfcv
CHILDREN
atomden
SOURCE
36 #if defined HAVE_CONFIG_H 37 #include "config.h" 38 #endif 39 40 #include "abi_common.h" 41 42 subroutine read_atomden(MPI_enreg,natom,nspden,ntypat,pawfgr, & 43 & rhor_paw,typat,rprimd,xred,prtvol,file_prefix) 44 45 use defs_basis 46 use defs_abitypes 47 use m_profiling_abi 48 use m_errors 49 50 use m_io_tools, only : open_file 51 use m_pawfgr, only : pawfgr_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 'read_atomden' 57 use interfaces_65_paw, except_this_one => read_atomden 58 !End of the abilint section 59 60 implicit none 61 62 !Arguments ------------------------------------ 63 !scalars 64 integer,intent(in) :: natom,nspden,ntypat,prtvol 65 type(pawfgr_type),intent(in) :: pawfgr 66 !arrays 67 type(MPI_type),intent(in) :: MPI_enreg 68 integer,intent(in) :: typat(natom) 69 real(dp), intent(in) :: rprimd(3,3),xred(3,natom) 70 real(dp),intent(inout) :: rhor_paw(pawfgr%nfft,nspden) 71 character(len=7), intent(in) :: file_prefix 72 73 !Local variables------------------------------- 74 !scalars 75 character(len=500) :: message 76 character(len=120) :: filename 77 character(len=7) :: calctype='replace' 78 integer :: igrid,i,i1,i2,i3,io_err,itypat,unt 79 integer :: natomgrmax,nlines,ngrid,n1,n2,n3 80 real(dp) :: difx,dify,difz,ucvol!,norm 81 !arrays 82 integer :: natomgr(ntypat) 83 real(dp) :: a(3),b(3),c(3) 84 real(dp) :: gmet(3,3),gprimd(3,3),rmet(3,3) 85 real(dp),allocatable :: atomrgrid(:,:),r_vec_grid(:,:),density(:,:) 86 real(dp),allocatable :: rho(:) 87 88 ! ************************************************************************ 89 90 !Initialise various variables 91 ngrid = pawfgr%nfft 92 a(:) = rprimd(:,1) 93 b(:) = rprimd(:,2) 94 c(:) = rprimd(:,3) 95 ABI_ALLOCATE(rho,(ngrid)) 96 if (nspden/=1) then 97 MSG_ERROR('read_atomden: Only nspden=1 allowed.') 98 end if 99 rho = rhor_paw(:,1) 100 gmet=zero;gprimd=zero;rmet=zero;ucvol=zero 101 102 103 !Calculate the r vector (reduced coord.) of the fine gridpoints 104 ABI_ALLOCATE(r_vec_grid,(3,ngrid)) 105 igrid = 0 106 n1 = pawfgr%ngfft(1) 107 n2 = pawfgr%ngfft(2) 108 n3 = pawfgr%ngfft(3) 109 do i3=0,n3-1 110 difz=dble(i3)/dble(n3) 111 do i2=0,n2-1 112 dify=dble(i2)/dble(n2) 113 do i1=0,n1-1 114 difx=dble(i1)/dble(n1) 115 igrid = igrid + 1 116 r_vec_grid(1,igrid)=difx*rprimd(1,1)+dify*rprimd(1,2)+difz*rprimd(1,3) 117 r_vec_grid(2,igrid)=difx*rprimd(2,1)+dify*rprimd(2,2)+difz*rprimd(2,3) 118 r_vec_grid(3,igrid)=difx*rprimd(3,1)+dify*rprimd(3,2)+difz*rprimd(3,3) 119 end do 120 end do 121 end do 122 if (igrid/=ngrid) then 123 MSG_ERROR('read_atomden: igrid not equal to ngrid') 124 end if 125 126 !Read in atomic density data for each atom type 127 !first check how many datapoints are in each file 128 do itypat=1,ntypat 129 filename='';io_err=0; 130 if (itypat>0) write(filename,'(a,a,i1,a)') trim(file_prefix), & 131 & '_density_atom_type',itypat,'.dat' 132 if (itypat>10) write(filename,'(a,a,i2,a)') trim(file_prefix), & 133 & '_density_atom_type',itypat,'.dat' 134 if (open_file(filename, message, newunit=unt, status='old',action='read') /= 0) then 135 write(std_out,*) 'ERROR in read_atomden: Could not open file: ',filename 136 write(std_out,*) ' Current implementation requires this file to be present' 137 write(std_out,*) ' for each type of atom.' 138 write(std_out,*)trim(message) 139 MSG_ERROR("Cannot continue") 140 end if 141 ! Check number of lines in file 142 nlines = 1;io_err=0; 143 do 144 read(unt,*,iostat=io_err) 145 if (io_err<0) exit 146 nlines = nlines + 1 147 end do 148 close(unt) 149 natomgr(itypat) = nlines - 2 150 end do ! Atom type 151 !Allocate arrays and read in data 152 natomgrmax = maxval(natomgr) 153 ABI_ALLOCATE(atomrgrid,(natomgrmax,ntypat)) 154 ABI_ALLOCATE(density,(natomgrmax,ntypat)) 155 atomrgrid = zero ; density = zero 156 do itypat=1,ntypat 157 filename='';io_err=0; 158 if (itypat>0) write(filename,'(a,a,i1,a)') trim(file_prefix), & 159 & '_density_atom_type',itypat,'.dat' 160 if (itypat>10) write(filename,'(a,a,i2,a)') trim(file_prefix), & 161 & '_density_atom_type',itypat,'.dat' 162 if (open_file(filename,message,newunit=unt,status='old',action='read') /= 0) then 163 MSG_ERROR(message) 164 end if 165 read(unt,*) ! Skip comment line 166 do i=1,natomgr(itypat) 167 read(unt,*) atomrgrid(i,itypat),density(i,itypat) 168 end do 169 close(unt) 170 if (atomrgrid(1,itypat)/=zero) then 171 write(std_out,*) 'ERROR in read_atomden, in file: ',filename 172 write(std_out,*) ' First gridpoint has to be the origin.' 173 MSG_ERROR("Cannot continue") 174 end if 175 end do ! Atom type 176 177 !write(std_out,*) '*** --- In read_atomden before call--- ***' 178 !write(std_out,*) ' calctype:',calctype,' natom:',natom 179 !write(std_out,*) ' ntypat:',ntypat,' typat:',typat 180 !write(std_out,*) ' ngrid:',ngrid 181 !write(std_out,*) ' a:',a 182 !write(std_out,*) ' b:',b 183 !write(std_out,*) ' c:',c 184 !write(std_out,*) ' xred:',xred 185 !write(std_out,*) ' natomgr:',natomgr 186 !write(std_out,*) 'natomgrmax:',natomgrmax 187 !write(std_out,*) ' atomrgrid:',atomrgrid 188 !write(std_out,*) ' density:',density 189 !write(std_out,*) 'r_vec_grid:' 190 !write(std_out,*) r_vec_grid 191 192 !Call atomden 193 call atomden(MPI_enreg,natom,ntypat,typat,ngrid,r_vec_grid,rho,a,b,c,xred, & 194 & natomgr,natomgrmax,atomrgrid,density,prtvol,calctype) 195 196 !if (prtvol>9) then ! calculate norm 197 !call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 198 !norm = SUM(rho(:))*ucvol/dble(n1*n2*n3) 199 !write(message,'(a,F8.4)') ' In read_atomden - NORM OF DENSITY: ',norm 200 !call wrtout(std_out,message,'COLL') 201 !end if 202 203 rhor_paw(:,1) = rho 204 205 if (allocated(atomrgrid)) then 206 ABI_DEALLOCATE(atomrgrid) 207 end if 208 if (allocated(density)) then 209 ABI_DEALLOCATE(density) 210 end if 211 if (allocated(r_vec_grid)) then 212 ABI_DEALLOCATE(r_vec_grid) 213 end if 214 if (allocated(rho)) then 215 ABI_DEALLOCATE(rho) 216 end if 217 218 return 219 220 end subroutine read_atomden