TABLE OF CONTENTS
ABINIT/symaxes [ Functions ]
NAME
symaxes
FUNCTION
Determines the type of symmetry operation, for the proper symmetries 2,2_1,3,3_1,3_2,4,4_1,4_2,4_3,6,6_1,...6_5
COPYRIGHT
Copyright (C) 2000-2018 ABINIT group (RC, 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
center=type of bravais lattice centering center=0 no centering center=-1 body-centered center=-3 face-centered center=1 A-face centered center=2 B-face centered center=3 C-face centered iholohedry=type of holohedry iholohedry=1 triclinic 1bar iholohedry=2 monoclinic 2/m iholohedry=3 orthorhombic mmm iholohedry=4 tetragonal 4/mmm iholohedry=5 trigonal 3bar m (rhombohedral Bravais latt) iholohedry=6 hexagonal 6/mmm iholohedry=7 cubic m3bar m isym=number of the symmetry operation that is currently analyzed isymrelconv=symrel matrix for the particular operation, in conv. axes ordersym=order of the symmetry operation tnons_order=order of the screw translation trialt(3)=screw translation associated with the symmetry operation in conventional axes (all components in the range ]-1/2,1/2] )
OUTPUT
label=a user friendly label for the rotation type_axis=type of the symmetry operation
NOTES
It is assumed that the symmetry operations will be entered in the symrel tnonsconv arrays, for the CONVENTIONAL cell. For proper symmetries (rotations), the associated translation is determined. There is a subtlety with translations associated with rotations : all the rotations with axis parallel to the one analysed do not all have the same translation characteristics. This is clearly seen in the extended Hermann-Mauguin symbols, see the international table for crystallography, chapter 4. In the treatment that we adopt, one will distinguish the cases of primitive Bravais lattices, and centered bravais lattices. In the latter case, in the present routine, at the exception of the trigonal axis for the cubic system, we explicitely generate the correct ratio of different translations, so that their type can be explicitely assigned, without confusion. By contrast, for primitive lattices, the "tnons" that has been transmitted to the present routine might be one of the few possible translations vectors, nearly at random. We deal with this case by the explicit examination of the system classes, and the identification of such a possibility. In particular: (1) for the trigonal axis in the rhombohedral Bravais lattice, or in the cubic system, there is an equal number of 3, 3_1, and 3_2 axes parallel to each other, in a cell that is primitive (as well as conventional). In this particular case, in the present routine, all 3, 3_1 and 3_2 axes are assigned to be 3 axes, independently of the centering. (2) for the 4- or 6- axes, no confusion is possible : in the primitive cell, there is only one possible translation, while in the centered cells, the correct ratio of translation vectors will be generated (3) for the binary axes, there is no problem when the cell is centered, but there are problems (3a) for the tP Bravais lattice, for an axis in a tertiary direction, (see the description of the lattice symmetry directions table 2.4.1 of the international tables for crystallography), where the family of axes is made equally of 2 and 2_1 axis. In this case, we attribute the binary axis to the specific class of "tertiary 2-axis". We keep track of the 2 or 2_1 characteristics of all other binary axes (3b) for the tI Bravais lattice, in all the directions, there is an equal number of 2 and 2_1 axes. We distinguish the primary and secondary family from the tertiary family. (3c) for the hP Bravais lattice, each binary axis can present no translation or be a screw axis (in the same direction). For primary axes, one need the "2" and "2_1" classification, while for secondary and tertiary axes, the associated translation vector will have not importance. However, one will need to distinguish secondary from tertiary, and these from primary axes. So, this is the most complicated case, for binary axes, with the following sets of binary axes : "2", "2_1", "secondary 2" and "tertiary 2". (3d) for the hR Bravais lattice, each binary axis can present no translation or be a screw axis (in the same direction). There is no distinction between tertiary axes and other, so that we simply assign a binary axis to "2-axis" (3e) for the cP lattice, the binary axes along tertiary directions can also have different translation vectors, while for the primary direction, there is no such ambiguity. So, we will attribute tertiary 2 axis to the "tertiary 2-axis" set (there are always 6), and attribute 2 and 2_1 primary axes to the corresponding sets.
PARENTS
symcharac
CHILDREN
wrtout
SOURCE
119 #if defined HAVE_CONFIG_H 120 #include "config.h" 121 #endif 122 123 #include "abi_common.h" 124 125 126 subroutine symaxes(center,iholohedry,isym,isymrelconv,label,ordersym,tnons_order,trialt,type_axis) 127 128 use defs_basis 129 use m_profiling_abi 130 use m_errors 131 132 !This section has been created automatically by the script Abilint (TD). 133 !Do not modify the following lines by hand. 134 #undef ABI_FUNC 135 #define ABI_FUNC 'symaxes' 136 use interfaces_14_hidewrite 137 !End of the abilint section 138 139 implicit none 140 141 !Arguments ------------------------------------ 142 !scalars 143 integer,intent(in) :: center,iholohedry,isym,ordersym,tnons_order 144 integer,intent(out) :: type_axis 145 character(len=128),intent(out) :: label 146 !arrays 147 integer,intent(in) :: isymrelconv(3,3) 148 real(dp),intent(in) :: trialt(3) 149 150 !Local variables------------------------------- 151 !scalars 152 logical,parameter :: verbose=.FALSE. 153 character(len=500) :: message 154 integer :: direction,directiontype 155 real(dp),parameter :: nzero=1.0d-6 156 157 !************************************************************************** 158 159 !DEBUG 160 !write(std_out,*)' symaxes : enter, isym=',isym 161 !write(std_out,*)' symaxes : iholohedry, ',iholohedry 162 !write(std_out,*)' symaxes : center, ',center 163 !stop 164 !ENDDEBUG 165 166 select case(ordersym) 167 168 case(2) ! point symmetry 2 169 ! Must characterize directiontype for cP, tP, tI, and hP Bravais lattices 170 directiontype=1 171 if( iholohedry==4 .or. iholohedry==7) then ! tP or cP Bravais lattices 172 if(abs(isymrelconv(1,1))+ & 173 & abs(isymrelconv(2,2))+ & 174 & abs(isymrelconv(3,3)) ==1) directiontype=3 175 else if(iholohedry==6)then ! hP Bravais lattice 176 if(sum(isymrelconv(:,:))/=-1 )directiontype=2 177 if(sum(isymrelconv(:,:))==0 .or. sum(isymrelconv(:,:))==-3 )& 178 & directiontype=3 179 ! directiontype=1 corresponds to a primary axis 180 ! directiontype=2 corresponds to a tertiary axis 181 ! directiontype=3 corresponds to a secondary axis 182 end if 183 184 ! DEBUG 185 ! write(std_out,*)' directiontype=',directiontype 186 ! write(std_out,'(a,3i6)' )' isymrelconv(1:3)=',isymrelconv(:,1) 187 ! write(std_out,'(a,3i6)' )' isymrelconv(4:6)=',isymrelconv(:,2) 188 ! write(std_out,'(a,3i6)' )' isymrelconv(7:9)=',isymrelconv(:,3) 189 ! write(std_out,'(a,i)' )' tnons_order=',tnons_order 190 ! ENDDEBUG 191 192 ! Now, classify the 2 axes 193 if(directiontype==2)then 194 type_axis=4 ! secondary 2 (only in the hP Bravais latt case) 195 write(label,'(a)') 'a secondary 2-axis ' 196 197 else if(directiontype==3 .and. iholohedry==4)then 198 type_axis=21 ! tertiary 2 199 write(label,'(a)') 'a tertiary 2-axis ' 200 else if(directiontype==3 .and. & 201 & center==0 .and. (iholohedry==6.or.iholohedry==7) )then 202 type_axis=21 ! tertiary 2 203 write(label,'(a)') 'a tertiary 2-axis ' 204 else if(tnons_order==1 .or. (iholohedry==4 .and. center==-1) .or. & 205 & iholohedry==5)then 206 type_axis=9 ! 2 207 write(label,'(a)') 'a 2-axis ' 208 else 209 type_axis=20 ! 2_1 210 write(label,'(a)') 'a 2_1-axis ' 211 end if 212 213 case(3) ! point symmetry 3 214 if(tnons_order==1)then 215 type_axis=10 ! 3 216 write(label,'(a)') 'a 3-axis ' 217 else if(iholohedry==5 .or. iholohedry==7)then 218 ! This is a special situation : in the same family of parallel 3-axis, 219 ! one will have an equal number of 3, 3_1 and 3_2 axes, so that 220 ! it is non-sense to try to classify one of them. 221 type_axis=10 ! 3, 3_1 or 3_2, undistinguishable 222 write(label,'(a)') 'a 3, 3_1 or 3_2 axis ' 223 else 224 ! DEBUG 225 ! write(std_out,*)'isymrelconv=',isymrelconv(:,:) 226 ! write(std_out,*)'trialt=',trialt(:) 227 ! ENDDEBUG 228 ! Must recognize 3_1 or 3_2 229 if(isymrelconv(1,1)==0)then ! 3+ 230 if(abs(trialt(3)-third)<nzero)type_axis=22 ! 3_1 231 if(abs(trialt(3)+third)<nzero)type_axis=23 ! 3_2 232 else if(isymrelconv(1,1)==-1)then ! 3- 233 if(abs(trialt(3)-third)<nzero)type_axis=23 ! 3_2 234 if(abs(trialt(3)+third)<nzero)type_axis=22 ! 3_1 235 end if 236 write(label,'(a)') 'a 3_1 or 3_2-axis ' 237 end if 238 239 case(4) ! point symmetry 4 240 if(tnons_order==1)then 241 type_axis=12 ! 4 242 write(label,'(a)') 'a 4-axis ' 243 else if(tnons_order==2)then 244 type_axis=25 ! 4_2 245 write(label,'(a)') 'a 4_2-axis ' 246 else if(center/=0)then 247 type_axis=24 ! 4_1 or 4_3 248 write(label,'(a)') 'a 4_1 or 4_3-axis ' 249 else 250 ! DEBUG 251 ! write(std_out,*)'isymrelconv=',isymrelconv(:,:) 252 ! write(std_out,*)'trialt=',trialt(:) 253 ! ENDDEBUG 254 ! Must recognize 4_1 or 4_3, along the three primary directions 255 do direction=1,3 256 if(isymrelconv(direction,direction)==1)then ! 257 if( (direction==1 .and. isymrelconv(2,3)==-1) .or. & 258 & (direction==2 .and. isymrelconv(3,1)==-1) .or. & 259 & (direction==3 .and. isymrelconv(1,2)==-1) )then ! 4+ 260 if(abs(trialt(direction)-quarter)<nzero)type_axis=24 ! 4_1 261 if(abs(trialt(direction)+quarter)<nzero)type_axis=26 ! 4_3 262 else if( (direction==1 .and. isymrelconv(2,3)==1) .or. & 263 & (direction==2 .and. isymrelconv(3,1)==1) .or. & 264 & (direction==3 .and. isymrelconv(1,2)==1) )then ! 4- 265 if(abs(trialt(direction)-quarter)<nzero)type_axis=26 ! 4_3 266 if(abs(trialt(direction)+quarter)<nzero)type_axis=24 ! 4_1 267 end if 268 end if 269 end do 270 write(label,'(a)') 'a 4_1 or 4_3-axis ' 271 end if 272 273 case(6) ! point symmetry 6 274 if(tnons_order==1)then 275 type_axis=14 ! 6 276 write(label,'(a)') 'a 6-axis ' 277 else if(tnons_order==2)then 278 type_axis=29 ! 6_3 279 write(label,'(a)') 'a 6_3-axis ' 280 else if(tnons_order==3)then 281 ! DEBUG 282 ! write(std_out,*)'isymrelconv=',isymrelconv(:,:) 283 ! write(std_out,*)'trialt=',trialt(:) 284 ! ENDDEBUG 285 ! Must recognize 6_2 or 6_4 286 if(isymrelconv(1,1)==1)then ! 6+ 287 if(abs(trialt(3)-third)<nzero)type_axis=28 ! 6_2 288 if(abs(trialt(3)+third)<nzero)type_axis=30 ! 6_4 289 else if(isymrelconv(1,1)==0)then ! 6- 290 if(abs(trialt(3)-third)<nzero)type_axis=30 ! 6_4 291 if(abs(trialt(3)+third)<nzero)type_axis=28 ! 6_2 292 end if 293 write(label,'(a)') 'a 6_2 or 6_4-axis ' 294 else 295 ! DEBUG 296 ! write(std_out,*)'isymrelconv=',isymrelconv(:,:) 297 ! write(std_out,*)'trialt=',trialt(:) 298 ! ENDDEBUG 299 ! Must recognize 6_1 or 6_5 300 if(isymrelconv(1,1)==1)then ! 6+ 301 if(abs(trialt(3)-sixth)<nzero)type_axis=27 ! 6_1 302 if(abs(trialt(3)+sixth)<nzero)type_axis=31 ! 6_5 303 else if(isymrelconv(1,1)==0)then ! 6- 304 if(abs(trialt(3)-sixth)<nzero)type_axis=31 ! 6_5 305 if(abs(trialt(3)+sixth)<nzero)type_axis=27 ! 6_1 306 end if 307 write(label,'(a)') 'a 6_1 or 6_5-axis ' 308 end if 309 310 end select 311 312 if (verbose) then 313 write(message,'(a,i3,a,a)')' symaxes : the symmetry operation no. ',isym,' is ', trim(label) 314 call wrtout(std_out,message,'COLL') 315 end if 316 317 end subroutine symaxes