TABLE OF CONTENTS
ABINIT/symbrav [ Functions ]
NAME
symbrav
FUNCTION
From the list of symmetry operations, and the lattice vectors, determine the Bravais information (including the holohedry, the centering, the coordinate of the primitive vectors in the conventional vectors), as well as the point group.
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
msym=dimension of symrel nsym=actual number of symmetries rprimd(3,3)=dimensional primitive translations for real space (bohr) symrel(3,3,msym)=symmetry operations in real space in terms of primitive translations tolsym=tolerance for the symmetries
OUTPUT
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)
ptgroup=symmetry point group
[axis(3)]=Invariant axis in the conventional vector coordinates
Set to (/0,0,0/) if the lattice belongs to the same holohedry as the lattice+atoms (+electric field + ...).
PARENTS
m_esymm,symanal
CHILDREN
matr3inv,symlatt,symptgroup,symrelrot
SOURCE
44 #if defined HAVE_CONFIG_H 45 #include "config.h" 46 #endif 47 48 #include "abi_common.h" 49 50 subroutine symbrav(bravais,msym,nsym,ptgroup,rprimd,symrel,tolsym,axis) 51 52 use defs_basis 53 use m_profiling_abi 54 use m_errors 55 56 !This section has been created automatically by the script Abilint (TD). 57 !Do not modify the following lines by hand. 58 #undef ABI_FUNC 59 #define ABI_FUNC 'symbrav' 60 use interfaces_32_util 61 use interfaces_41_geometry, except_this_one => symbrav 62 !End of the abilint section 63 64 implicit none 65 66 !Arguments ------------------------------------ 67 !scalars 68 integer,intent(in) :: msym,nsym 69 real(dp),intent(in) :: tolsym 70 character(len=5),intent(out) :: ptgroup 71 !arrays 72 integer,intent(in) :: symrel(3,3,msym) 73 integer,optional,intent(out) :: axis(3) 74 integer,intent(out) :: bravais(11) 75 real(dp),intent(in) :: rprimd(3,3) 76 77 !Local variables------------------------------- 78 !scalars 79 integer :: iaxis,ii,bravais1now,ideform,iholohedry,invariant,isym 80 integer :: jaxis,next_stage,nptsym,problem 81 real(dp) :: norm,scprod 82 character(len=500) :: message 83 !arrays 84 integer :: identity(3,3),axis_trial(3),hexa_axes(3,7),ortho_axes(3,13) 85 integer,allocatable :: ptsymrel(:,:,:),symrelconv(:,:,:) 86 real(dp) :: axes(3,3),axis_cart(3),axis_red(3) 87 real(dp) :: rprimdconv(3,3),rprimdtry(3,3),rprimdnow(3,3) 88 real(dp) :: rprimdconv_invt(3,3) 89 90 !************************************************************************** 91 92 identity(:,:)=0 93 identity(1,1)=1 ; identity(2,2)=1 ; identity(3,3)=1 94 95 ortho_axes(:,:)=0 96 ortho_axes(1,1)=1 97 ortho_axes(2,2)=1 98 ortho_axes(3,3)=1 99 ortho_axes(:,4)=(/0,1,1/) 100 ortho_axes(:,5)=(/1,0,1/) 101 ortho_axes(:,6)=(/1,1,0/) 102 ortho_axes(:,7)=(/0,1,-1/) 103 ortho_axes(:,8)=(/-1,0,1/) 104 ortho_axes(:,9)=(/1,-1,0/) 105 ortho_axes(:,10)=(/1,1,1/) 106 ortho_axes(:,11)=(/-1,1,1/) 107 ortho_axes(:,12)=(/1,-1,1/) 108 ortho_axes(:,13)=(/1,1,-1/) 109 110 hexa_axes(:,:)=0 111 hexa_axes(1,1)=1 112 hexa_axes(2,2)=1 113 hexa_axes(3,3)=1 114 hexa_axes(:,4)=(/1,-1,0/) 115 hexa_axes(:,5)=(/2,1,0/) 116 hexa_axes(:,6)=(/1,1,0/) 117 hexa_axes(:,7)=(/1,2,0/) 118 119 !Determine the point group from the list of symmetry operations. 120 !Also determine the holohedry, up to one undeterminacy : hR versus hP 121 call symptgroup(iholohedry,nsym,ptgroup,symrel) 122 123 !Loop over trial deformations 124 !This is needed in case the Bravais lattice determination from the lattice vectors 125 !has a higher holohedry than the real one, in which the symmetry 126 !operations for the atoms (or electric field, etc) are taken into account 127 iaxis=0 128 invariant=0 129 next_stage=0 130 rprimdnow(:,:)=rprimd(:,:) 131 rprimdtry(:,:)=rprimd(:,:) 132 ABI_ALLOCATE(symrelconv,(3,3,nsym)) 133 134 !At most will have to try 65 deformations (13 axes, five stages) 135 do ideform=1,65 136 137 ABI_ALLOCATE(ptsymrel,(3,3,msym)) 138 call symlatt(bravais,msym,nptsym,ptsymrel,rprimdtry,tolsym) 139 ABI_DEALLOCATE(ptsymrel) 140 141 ! Examine the agreement with bravais(1) 142 ! Warning : might change Bravais lattice hR to hP, if hexagonal axes 143 problem=0 144 select case (bravais(1)) 145 case(7) 146 if(iholohedry<6)problem=1 147 if(iholohedry==6)problem=2 148 case(6) 149 if(iholohedry<4)problem=1 150 if(iholohedry==7 .or. iholohedry==4)problem=2 151 ! Here, change hR into hP 152 if(iholohedry==5)iholohedry=6 153 case(5) 154 if(iholohedry<4)problem=1 155 if(iholohedry==7 .or. iholohedry==6 .or. iholohedry==4)problem=2 156 case(4) 157 if(iholohedry<4)problem=1 158 if(iholohedry>4)problem=2 159 case(3) 160 if(iholohedry<3)problem=1 161 if(iholohedry>3)problem=2 162 case(2) 163 if(iholohedry<2)problem=1 164 if(iholohedry>2)problem=2 165 case(1) 166 if(iholohedry>1)problem=2 167 end select 168 169 ! This is the usual situation, in which the lattice belong to the same holohedry 170 ! as the lattice+atoms (+electric field + ...) 171 if(problem==0)exit 172 173 if(problem==2)then 174 if(iaxis==0)then 175 write(message, '(3a,i3,3a,i3,7a)' )& 176 & 'The Bravais lattice determined only from the primitive',ch10,& 177 & 'vectors (rprim or angdeg), bravais(1)=',bravais(1),', is not compatible',ch10,& 178 & 'with the real one, iholohedry=',iholohedry,', obtained by taking into',ch10,& 179 & 'account the symmetry operations. This might be due to an insufficient',ch10,& 180 & 'number of digits in the specification of rprim (at least 10),',ch10,& 181 & 'or to an erroneous rprim or angdeg. If this is not the case, then ...' 182 MSG_BUG(message) 183 end if 184 if(iaxis==1)then 185 write(message, '(3a,3i3,2a,i3,2a,i3)' )& 186 & 'Could not succeed to determine the bravais lattice',ch10,& 187 & 'problem,iaxis,invariant=',problem,iaxis,invariant,ch10,& 188 & 'bravais(1)=',bravais(1),ch10,& 189 & 'iholohedry=',iholohedry 190 MSG_BUG(message) 191 end if 192 end if 193 194 if(problem==1)then ! One is left with the problem=1 case, basically iholohedry is lower than bravais(1) 195 if(iaxis==0)then 196 write(message, '(a,a,a,i3,a,a,a,i3,a,a,a)' )& 197 & 'The Bravais lattice determined only from the primitive',ch10,& 198 & 'vectors, bravais(1)=',bravais(1),', is more symmetric',ch10,& 199 & 'than the real one, iholohedry=',iholohedry,', obtained by taking into',ch10,& 200 & 'account the atomic positions. Start deforming the primitive vector set.' 201 MSG_COMMENT(message) 202 next_stage=1 203 else if(iaxis/=0)then 204 if(bravais(1)<bravais1now)then 205 write(message, '(3a,i3,3a,i3,2a)' )& 206 & 'The Bravais lattice determined from modified primitive',ch10,& 207 & 'vectors, bravais(1)=',bravais(1),', has a lower symmetry than before,',ch10,& 208 & 'but is still more symmetric than the real one, iholohedry=',iholohedry,ch10,& 209 & 'obtained by taking into account the atomic positions.' 210 MSG_COMMENT(message) 211 next_stage=1 212 else if(iaxis==1)then 213 write(message, '(3a,3i3,2a,i3,2a,i3)' )& 214 & 'Could not succeed to determine the bravais lattice',ch10,& 215 & 'problem,iaxis,invariant=',problem,iaxis,invariant,ch10,& 216 & 'bravais(1)=',bravais(1),ch10,& 217 & 'iholohedry=',iholohedry 218 MSG_BUG(message) 219 end if 220 end if 221 end if ! problem==1 222 223 if(next_stage==1)then 224 bravais1now=bravais(1) 225 rprimdnow(:,:)=rprimdtry(:,:) 226 ! Generate the symmetry operations in the conventional vector coordinates 227 rprimdconv(:,1)=bravais(3:5) 228 rprimdconv(:,2)=bravais(6:8) 229 rprimdconv(:,3)=bravais(9:11) 230 axes(:,:)=zero 231 axes(1,1)=one ; axes(2,2)=one ; axes(3,3)=one 232 symrelconv(:,:,1:nsym)=symrel(:,:,1:nsym) 233 call symrelrot(nsym,rprimdconv,axes,symrelconv,tolsym) 234 if(bravais(1)/=6)then 235 iaxis=14 236 else 237 iaxis=8 238 end if 239 next_stage=0 240 end if 241 242 iaxis=iaxis-1 243 do jaxis=iaxis,1,-1 244 if(bravais(1)/=6)then 245 axis_trial(:)=ortho_axes(:,jaxis) 246 else 247 axis_trial(:)=hexa_axes(:,jaxis) 248 end if 249 ! DEBUG 250 ! write(std_out,*)' symbrav : try jaxis=',jaxis 251 ! write(std_out,*)' axis_trial=',axis_trial 252 ! ENDDEBUG 253 invariant=1 254 ! Examine whether all symmetry operations leave the axis invariant (might be reversed, though) 255 do isym=1,nsym 256 if(sum(abs(matmul(symrelconv(:,:,isym),axis_trial)+(-axis_trial(:))))/=0 .and. & 257 & sum(abs(matmul(symrelconv(:,:,isym),axis_trial)+axis_trial(:)))/=0 )invariant=0 258 end do 259 if(invariant==1)then 260 iaxis=jaxis 261 ! write(message, '(2a,i3)' )ch10,' symbrav : found invariant axis, jaxis=',iaxis 262 ! call wrtout(std_out,message,'COLL') 263 exit 264 end if 265 end do 266 267 if(invariant==0)then 268 ! Not a single axis was invariant with respect to all operations ?! 269 ! do isym=1,nsym; write(std_out, '(a,10i4)' )' isym,symrelconv=',isym,symrelconv(:,:,isym); enddo 270 write(message, '(3a,3i3,2a,i3,2a,i3)' )& 271 & 'Could not succeed to determine the bravais lattice (not a single invariant)',ch10,& 272 & 'problem,iaxis,invariant=',problem,iaxis,invariant,ch10,& 273 & 'bravais(1)=',bravais(1),ch10,& 274 & 'iholohedry=',iholohedry 275 MSG_BUG(message) 276 end if 277 278 call matr3inv(rprimdconv,rprimdconv_invt) 279 axis_red(:)=axis_trial(1)*rprimdconv_invt(1,:)+ & 280 & axis_trial(2)*rprimdconv_invt(2,:)+ & 281 & axis_trial(3)*rprimdconv_invt(3,:) 282 axis_cart(:)=axis_red(1)*rprimdnow(:,1)+ & 283 & axis_red(2)*rprimdnow(:,2)+ & 284 & axis_red(3)*rprimdnow(:,3) 285 norm=sum(axis_cart(:)**2) 286 ! Expand by a uniform, quite arbitrary, dilatation, along the invariant axis 287 ! Note : make these dilatation different, according to ideform 288 ! XG 20151221 : Still, the interplay between the size of the deformation and the tolsym is not easy to address. 289 ! Indeed the deformation must be sufficiently large to be perceived by symlatt as a real breaking of the 290 ! symmetry of the lattice. In order to deal with all the small values od tolsym, it has been set at a minimum of tol3, 291 ! but it must also be larger than tolsym. Moreover, for some axis choice, the deformation is not aligned with the axis, decreasing 292 ! the effective deformation length. An additional factor of three is thus included, actually increased to six just to be sure... 293 do ii=1,3 294 scprod=axis_cart(1)*rprimdnow(1,ii)+axis_cart(2)*rprimdnow(2,ii)+axis_cart(3)*rprimdnow(3,ii) 295 rprimdtry(:,ii)=rprimdnow(:,ii)+ideform*(max(tol3,six*tolsym)-tol6)*scprod/norm*axis_cart(:) 296 end do 297 298 end do ! ideform 299 300 if(bravais(1)/=iholohedry)then 301 write(message, '(3a,3i3,2a,i3,2a,i3)' )& 302 & 'Despite efforts, Could not succeed to determine the bravais lattice :',ch10,& 303 & 'bravais(1)=',bravais(1),ch10,& 304 & 'iholohedry=',iholohedry 305 MSG_BUG(message) 306 end if 307 308 ABI_DEALLOCATE(symrelconv) 309 310 if (PRESENT(axis)) then ! Return symmetry axis. 311 axis=(/0,0,0/) 312 if (iaxis/=0) then 313 if(bravais(1)/=6)then 314 axis=ortho_axes(:,iaxis) 315 else 316 axis=hexa_axes(:,iaxis) 317 end if 318 end if 319 end if 320 321 end subroutine symbrav