TABLE OF CONTENTS
ABINIT/invacuum [ Functions ]
NAME
invacuum
FUNCTION
Determine whether there is vacuum along some of the primitive directions in real space.
COPYRIGHT
Copyright (C) 2000-2018 ABINIT group (XG) 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
jdtset=number of the dataset looked for lenstr=actual length of the string natom=number of atoms rprimd(3,3)=dimensional real space primitive translations (bohr) string*(*)=character string containing all the input data. Initialized previously in instrng. xred(3,natom)=reduced dimensionless atomic coordinates
OUTPUT
vacuum(3)= for each direction, 0 if no vacuum, 1 if vacuum
PARENTS
invars1,invars2
CHILDREN
intagm,metric,sort_dp
SOURCE
38 #if defined HAVE_CONFIG_H 39 #include "config.h" 40 #endif 41 42 #include "abi_common.h" 43 44 subroutine invacuum(jdtset,lenstr,natom,rprimd,string,vacuum,xred) 45 46 use defs_basis 47 use m_profiling_abi 48 use m_sort 49 50 use m_geometry, only : metric 51 use m_parser, only : intagm 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 'invacuum' 57 !End of the abilint section 58 59 implicit none 60 61 !Arguments ------------------------------------ 62 !scalars 63 integer,intent(in) :: jdtset,lenstr,natom 64 character(len=*),intent(in) :: string 65 !arrays 66 integer,intent(out) :: vacuum(3) 67 real(dp),intent(in) :: rprimd(3,3),xred(3,natom) 68 69 !Local variables------------------------------- 70 !scalars 71 integer :: ia,ii,marr,tread 72 real(dp) :: max_diff_xred,ucvol,vacwidth,vacxred 73 !arrays 74 integer,allocatable :: list(:) 75 integer,allocatable :: intarr(:) 76 real(dp) :: gmet(3,3),gprimd(3,3),rmet(3,3) 77 real(dp),allocatable :: xred_sorted(:) 78 real(dp),allocatable :: dprarr(:) 79 80 ! ************************************************************************* 81 82 !Compute the maximum size of arrays intarr and dprarr 83 marr=3 84 ABI_ALLOCATE(intarr,(marr)) 85 ABI_ALLOCATE(dprarr,(marr)) 86 87 !Get metric quantities 88 call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 89 90 !Read vacwidth, or set the default 91 vacwidth=10.0_dp 92 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vacwidth',tread,'LEN') 93 if(tread==1) vacwidth=dprarr(1) 94 95 !Read vacuum, or compute it using the atomic coordinates and vacwidth. 96 vacuum(1:3)=0 97 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'vacuum',tread,'INT') 98 99 if(tread==1)then 100 vacuum(1:3)=intarr(1:3) 101 else 102 ! For each direction, determine whether a vacuum space exists 103 ABI_ALLOCATE(list,(natom)) 104 ABI_ALLOCATE(xred_sorted,(natom)) 105 do ii=1,3 106 ! This is the minimum xred difference needed to have vacwidth 107 vacxred=vacwidth*sqrt(sum(gprimd(:,ii)**2)) 108 ! Project the reduced coordinate in the [0.0_dp,1.0_dp[ interval 109 xred_sorted(:)=mod(xred(ii,:),1.0_dp) 110 ! list is dummy 111 list(:)=0 112 ! Sort xred_sorted 113 call sort_dp(natom,xred_sorted,list,tol14) 114 if(natom==1)then 115 max_diff_xred=1.0_dp 116 else 117 ! Compute the difference between each pair of atom in the sorted order 118 max_diff_xred=0.0_dp 119 do ia=1,natom-1 120 max_diff_xred=max(max_diff_xred,xred_sorted(ia+1)-xred_sorted(ia)) 121 end do 122 ! Do not forget the image of the first atom in the next cell 123 max_diff_xred=max(max_diff_xred,1.0_dp+xred_sorted(1)-xred_sorted(ia)) 124 end if 125 if(vacxred<max_diff_xred+tol10)vacuum(ii)=1 126 end do 127 ABI_DEALLOCATE(list) 128 ABI_DEALLOCATE(xred_sorted) 129 end if 130 131 !DEBUG 132 !write(std_out,*)' invacuum : vacuum=',vacuum(1:3) 133 !ENDDEBUG 134 135 ABI_DEALLOCATE(intarr) 136 ABI_DEALLOCATE(dprarr) 137 138 end subroutine invacuum