TABLE OF CONTENTS
ABINIT/symmetrize_rprimd [ Functions ]
NAME
symmetrize_rprimd
FUNCTION
Supposing the input rprimd does not preserve the length and angles following the symmetries, will generates a new set rprimd, on the basis of the expected characteristics of the conventional cell, as specified in bravais(:)
COPYRIGHT
Copyright (C) 2015-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
bravais(11): bravais(1)=iholohedry bravais(2)=center bravais(3:11)=coordinates of rprimd in the axes of the conventional bravais lattice (*2 if center/=0) nsym=actual number of symmetries symrel(3,3,1:nsym)=symmetry operations in real space in terms of primitive translations tolsym=tolerance for the symmetry operations (only for checking purposes, the new set rprimd will be coherent with the symmetry operations at a much accurate level).
SIDE EFFECTS
rprimd(3,3)=dimensional primitive translations for real space (bohr)
PARENTS
ingeo
CHILDREN
chkorthsy,holocell,matr3inv,metric
SOURCE
40 #if defined HAVE_CONFIG_H 41 #include "config.h" 42 #endif 43 44 #include "abi_common.h" 45 46 47 subroutine symmetrize_rprimd(bravais,nsym,rprimd,symrel,tolsym) 48 49 use defs_basis 50 use m_errors 51 use m_profiling_abi 52 53 use m_geometry, only : metric 54 55 !This section has been created automatically by the script Abilint (TD). 56 !Do not modify the following lines by hand. 57 #undef ABI_FUNC 58 #define ABI_FUNC 'symmetrize_rprimd' 59 use interfaces_32_util 60 use interfaces_41_geometry, except_this_one => symmetrize_rprimd 61 !End of the abilint section 62 63 implicit none 64 65 !Arguments ------------------------------------ 66 !scalars 67 integer,intent(in) :: nsym 68 real(dp),intent(in) :: tolsym 69 !arrays 70 integer,intent(in) :: bravais(11),symrel(3,3,nsym) 71 real(dp),intent(inout) :: rprimd(3,3) 72 73 !Local variables------------------------------- 74 !scalars 75 integer :: foundc,iexit,ii,jj 76 real(dp):: ucvol,reldiff 77 character(len=500) :: msg 78 !arrays 79 real(dp):: aa(3,3),ait(3,3),cell_base(3,3),gmet(3,3),gprimd(3,3),rmet(3,3),rprimd_new(3,3) 80 81 ! ************************************************************************* 82 83 !Build the conventional cell basis vectors in cartesian coordinates 84 aa(:,1)=bravais(3:5) 85 aa(:,2)=bravais(6:8) 86 aa(:,3)=bravais(9:11) 87 !Inverse transpose 88 call matr3inv(aa,ait) 89 do ii=1,3 90 cell_base(:,ii)=ait(ii,1)*rprimd(:,1)+ait(ii,2)*rprimd(:,2)+ait(ii,3)*rprimd(:,3) 91 end do 92 93 !Enforce the proper holohedry on the conventional cell vectors. 94 call holocell(cell_base,1,foundc,bravais(1),tolsym) 95 96 !Reconstruct the dimensional primitive vectors 97 do ii=1,3 98 rprimd_new(:,ii)=aa(1,ii)*cell_base(:,1)+aa(2,ii)*cell_base(:,2)+aa(3,ii)*cell_base(:,3) 99 end do 100 101 !Check whether the modification make sense 102 do ii=1,3 103 do jj=1,3 104 reldiff=(rprimd_new(ii,jj)-rprimd(ii,jj))/sqrt(sum(rprimd(:,jj)**2)) 105 ! Allow for twice tolsym 106 if(abs(reldiff)>two*tolsym)then 107 write(msg,'(a,6(2a,3es14.6))')& 108 & 'Failed rectification of lattice vectors to comply with Bravais lattice identification, modifs are too large',ch10,& 109 & ' rprimd =',rprimd(:,1),ch10,& 110 & ' ',rprimd(:,2),ch10,& 111 & ' ',rprimd(:,3),ch10,& 112 & ' rprimd_new=',rprimd_new(:,1),ch10,& 113 & ' ',rprimd_new(:,2),ch10,& 114 & ' ',rprimd_new(:,3) 115 MSG_ERROR_CLASS(msg, "TolSymError") 116 end if 117 end do 118 end do 119 120 rprimd(:,:)=rprimd_new(:,:) 121 122 !Check whether the symmetry operations are consistent with the lattice vectors 123 call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 124 iexit=1 125 126 call chkorthsy(gprimd,iexit,nsym,rmet,rprimd,symrel) 127 128 end subroutine symmetrize_rprimd
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-2018 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
178 subroutine symmetrize_xred(indsym,natom,nsym,symrel,tnons,xred) 179 180 use defs_basis 181 use m_profiling_abi 182 183 !This section has been created automatically by the script Abilint (TD). 184 !Do not modify the following lines by hand. 185 #undef ABI_FUNC 186 #define ABI_FUNC 'symmetrize_xred' 187 !End of the abilint section 188 189 implicit none 190 191 !Arguments ------------------------------------ 192 !scalars 193 integer,intent(in) :: natom,nsym 194 !arrays 195 integer,intent(in) :: indsym(4,nsym,natom),symrel(3,3,nsym) 196 real(dp),intent(in) :: tnons(3,nsym) 197 real(dp),intent(inout) :: xred(3,natom) 198 199 !Local variables------------------------------- 200 !scalars 201 integer :: iatom,ib,isym 202 integer :: ii,jj 203 real(dp) :: fc1,fc2,fc3 204 real(dp) :: diff 205 logical :: dissimilar 206 !arrays 207 real(dp) :: tsum(3),tt(3) 208 real(dp),allocatable :: xredsym(:,:) 209 real(dp) :: transl(3) ! translation vector 210 211 ! ************************************************************************* 212 ! 213 !Check whether group contains more than identity; 214 !if not then simply return 215 if (nsym>1) then 216 217 ! loop over atoms 218 ABI_ALLOCATE(xredsym,(3,natom)) 219 do iatom=1,natom 220 tsum(1)=0.0d0 221 tsum(2)=0.0d0 222 tsum(3)=0.0d0 223 ! 224 ! loop over symmetries 225 do isym=1,nsym 226 ! atom ib is atom into which iatom is rotated by inverse of 227 ! symmetry isym (inverse of symrel(mu,nu,isym)) 228 ib=indsym(4,isym,iatom) 229 ! Find the reduced coordinates after translation=t(indsym)+transl 230 fc1=xred(1,ib)+dble(indsym(1,isym,iatom)) 231 fc2=xred(2,ib)+dble(indsym(2,isym,iatom)) 232 fc3=xred(3,ib)+dble(indsym(3,isym,iatom)) 233 ! Compute [S * (x(indsym)+transl) ] + tnonsymmorphic 234 tt(:)=dble(symrel(:,1,isym))*fc1+& 235 & dble(symrel(:,2,isym))*fc2+& 236 & dble(symrel(:,3,isym))*fc3+ tnons(:,isym) 237 238 ! Average over nominally equivalent atomic positions 239 tsum(:)=tsum(:)+tt(:) 240 end do 241 ! 242 ! Set symmetrized result to sum over number of terms 243 xredsym(:,iatom)=tsum(:)/dble(nsym) 244 245 ! End loop over iatom 246 end do 247 248 transl(:)=xredsym(:,1)-nint(xredsym(:,1)) 249 250 ! Compute the smallest translation to an integer 251 do jj=2,natom 252 do ii=1,3 253 diff=xredsym(ii,jj)-nint(xredsym(ii,jj)) 254 if (diff<transl(ii)) transl(ii)=diff 255 end do 256 end do 257 258 ! Test if the translation on each direction is small 259 ! Tolerance 1E-13 260 do ii=1,3 261 if (abs(transl(ii))>1e-13) transl(ii)=0.0 262 end do 263 264 ! Execute translation 265 do jj=1,natom 266 do ii=1,3 267 xredsym(ii,jj)=xredsym(ii,jj)-transl(ii) 268 end do 269 end do 270 271 ! Test if xredsym is too similar to xred 272 ! Tolerance 1E-15 273 dissimilar=.FALSE. 274 do jj=1,natom 275 do ii=1,3 276 if (abs(xredsym(ii,jj)-xred(ii,jj))>1E-15) dissimilar=.TRUE. 277 end do 278 end do 279 280 if (dissimilar) xred(:,:)=xredsym(:,:) 281 ABI_DEALLOCATE(xredsym) 282 283 ! End condition of nsym/=1 284 end if 285 286 end subroutine symmetrize_xred