TABLE OF CONTENTS
ABINIT/symmetrize_xred [ Functions ]
NAME
symmetrize_xred
FUNCTION
Symmetrize atomic coordinates using input symmetry matrices symrel which are expressed in terms of the basis of real space primitive translations (array elements are integers). Input array indsym(4,isym,iatom) gives label of atom into which iatom is rotated by INVERSE of symmetry element isym and also gives primitive translation to get back to unit cell. This version uses improvement in algorithm suggested by Andrew Horsfield (see symatm.f).
COPYRIGHT
Copyright (C) 1998-2017 ABINIT group (DCA, XG, GMR) 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
indsym(4,nsym,natom)=indirect indexing array giving label of atom into which iatom is rotated by symmetry element isym natom=number of atoms nsym=number of symmetries in group symrel(3,3,nsym)=symmetry matrices in terms of real space primitive translations tnons(3,nsym)=nonsymmorphic translations for symmetries
OUTPUT
(see side effects)
SIDE EFFECTS
Input/Output xred(3,natom)= (input) atomic coordinates in terms of real space translations (output) symmetrized atomic coordinates in terms of real space translations
PARENTS
ingeo,mover,nonlinear,respfn,scfcv
CHILDREN
SOURCE
49 #if defined HAVE_CONFIG_H 50 #include "config.h" 51 #endif 52 53 #include "abi_common.h" 54 55 subroutine symmetrize_xred(indsym,natom,nsym,symrel,tnons,xred) 56 57 use defs_basis 58 use m_profiling_abi 59 60 !This section has been created automatically by the script Abilint (TD). 61 !Do not modify the following lines by hand. 62 #undef ABI_FUNC 63 #define ABI_FUNC 'symmetrize_xred' 64 !End of the abilint section 65 66 implicit none 67 68 !Arguments ------------------------------------ 69 !scalars 70 integer,intent(in) :: natom,nsym 71 !arrays 72 integer,intent(in) :: indsym(4,nsym,natom),symrel(3,3,nsym) 73 real(dp),intent(in) :: tnons(3,nsym) 74 real(dp),intent(inout) :: xred(3,natom) 75 76 !Local variables------------------------------- 77 !scalars 78 integer :: iatom,ib,isym 79 integer :: ii,jj 80 real(dp) :: fc1,fc2,fc3 81 real(dp) :: diff 82 logical :: dissimilar 83 !arrays 84 real(dp) :: tsum(3),tt(3) 85 real(dp),allocatable :: xredsym(:,:) 86 real(dp) :: transl(3) ! translation vector 87 88 ! ************************************************************************* 89 ! 90 !Check whether group contains more than identity; 91 !if not then simply return 92 if (nsym>1) then 93 94 ! loop over atoms 95 ABI_ALLOCATE(xredsym,(3,natom)) 96 do iatom=1,natom 97 tsum(1)=0.0d0 98 tsum(2)=0.0d0 99 tsum(3)=0.0d0 100 ! 101 ! loop over symmetries 102 do isym=1,nsym 103 ! atom ib is atom into which iatom is rotated by inverse of 104 ! symmetry isym (inverse of symrel(mu,nu,isym)) 105 ib=indsym(4,isym,iatom) 106 ! Find the reduced coordinates after translation=t(indsym)+transl 107 fc1=xred(1,ib)+dble(indsym(1,isym,iatom)) 108 fc2=xred(2,ib)+dble(indsym(2,isym,iatom)) 109 fc3=xred(3,ib)+dble(indsym(3,isym,iatom)) 110 ! Compute [S * (x(indsym)+transl) ] + tnonsymmorphic 111 tt(:)=dble(symrel(:,1,isym))*fc1+& 112 & dble(symrel(:,2,isym))*fc2+& 113 & dble(symrel(:,3,isym))*fc3+ tnons(:,isym) 114 115 ! Average over nominally equivalent atomic positions 116 tsum(:)=tsum(:)+tt(:) 117 end do 118 ! 119 ! Set symmetrized result to sum over number of terms 120 xredsym(:,iatom)=tsum(:)/dble(nsym) 121 122 ! End loop over iatom 123 end do 124 125 transl(:)=xredsym(:,1)-nint(xredsym(:,1)) 126 127 ! Compute the smallest translation to an integer 128 do jj=2,natom 129 do ii=1,3 130 diff=xredsym(ii,jj)-nint(xredsym(ii,jj)) 131 if (diff<transl(ii)) transl(ii)=diff 132 end do 133 end do 134 135 ! Test if the translation on each direction is small 136 ! Tolerance 1E-13 137 do ii=1,3 138 if (abs(transl(ii))>1e-13) transl(ii)=0.0 139 end do 140 141 ! Execute translation 142 do jj=1,natom 143 do ii=1,3 144 xredsym(ii,jj)=xredsym(ii,jj)-transl(ii) 145 end do 146 end do 147 148 ! Test if xredsym is too similar to xred 149 ! Tolerance 1E-15 150 dissimilar=.FALSE. 151 do jj=1,natom 152 do ii=1,3 153 if (abs(xredsym(ii,jj)-xred(ii,jj))>1E-15) dissimilar=.TRUE. 154 end do 155 end do 156 157 if (dissimilar) xred(:,:)=xredsym(:,:) 158 ABI_DEALLOCATE(xredsym) 159 160 ! End condition of nsym/=1 161 end if 162 163 end subroutine symmetrize_xred