TABLE OF CONTENTS
ABINIT/symspgr [ Functions ]
NAME
symspgr
FUNCTION
Using the type of each symmetry operation (found in symplanes.f and symaxes.f): proper symmetries 1,2,2_1,3,3_1,3_2,4,4_1,4_2,4_3,6,6_1,...6_5 improper symmetries -1,m,a,b,c,d,n,g,-3,-4,-6 , build an array with the number of such operations. then, call symlist.f to identify the space group. The identification is not unambiguous still ...
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
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,nsym)= nsym symmetry operations in real space in terms of primitive translations tnons(3,nsym)=nonsymmorphic translations for each symmetry (would be 0 0 0 each for a symmorphic space group)
OUTPUT
spgroup=symmetry space group number
NOTES
It is assumed that the symmetry operations will be entered in the symrel tnons arrays, for the PRIMITIVE cell. The matrix of transformation from the primitive cell to the conventional cell is described in the array "bravais" (see symlatt.F90). The present routine first make the transformation from the primitive coordinates to the conventional ones, then eventually generate additional symmetries, taking into account the centering translations. Then, the order and determinant of each symmetry operation is determined. For proper symmetries (rotations), the associated translation is also determined. However, left or right handed screw rotations are not (presently) distinguished, and will be attributed equally to left or right. For the detailed description of the labelling of the axes, see symaxes.f and symplanes.f
PARENTS
symanal
CHILDREN
spgdata,symcharac,symdet,symlist_bcc,symlist_fcc,symlist_others symlist_prim,symrelrot,wrtout,xred2xcart
SOURCE
64 #if defined HAVE_CONFIG_H 65 #include "config.h" 66 #endif 67 68 #include "abi_common.h" 69 70 subroutine symspgr(bravais,nsym,spgroup,symrel,tnons,tolsym) 71 72 use defs_basis 73 use m_errors 74 use m_profiling_abi 75 use m_numeric_tools, only : OPERATOR(.x.) 76 77 use m_geometry, only : xred2xcart 78 79 !This section has been created automatically by the script Abilint (TD). 80 !Do not modify the following lines by hand. 81 #undef ABI_FUNC 82 #define ABI_FUNC 'symspgr' 83 use interfaces_14_hidewrite 84 use interfaces_41_geometry, except_this_one => symspgr 85 !End of the abilint section 86 87 implicit none 88 89 !Arguments ------------------------------------ 90 !scalars 91 integer,intent(in) :: nsym 92 integer,intent(out) :: spgroup 93 real(dp),intent(in) :: tolsym 94 !arrays 95 integer,intent(in) :: bravais(11),symrel(3,3,nsym) 96 real(dp),intent(inout) :: tnons(3,nsym) 97 98 !Local variables------------------------------- 99 !scalars 100 ! logical,parameter :: verbose=.FALSE. 101 integer :: additional_info,brvltt,center,direction=0,found,iholohedry,ii 102 integer :: ishift,isym,jj,nshift,nsymconv,spgaxor,spgorig,sporder 103 character(len=1) :: brvsb 104 character(len=15) :: intsb,ptintsb,ptschsb,schsb 105 character(len=35) :: intsbl 106 character(len=500) :: message 107 character(len=128) :: label 108 !arrays 109 integer :: n_axes(31),n_axest(31),prime(5),test_direction(3),symrel_uni(3,3) 110 integer :: uniaxis(3),uniaxis_try(3) 111 integer,allocatable :: determinant(:),symrelconv(:,:,:),t_axes(:) 112 real(dp) :: axes(3,3),rprimdconv(3,3),trialt(3),vect(3,3) 113 real(dp),allocatable :: shift(:,:),tnonsconv(:,:) 114 115 !************************************************************************** 116 117 DBG_ENTER("COLL") 118 119 !Initialize brvltt, from bravais(2) and bravais(1) 120 center=bravais(2) 121 iholohedry=bravais(1) 122 brvltt=1 123 if(center==-1)brvltt=2 ! Inner centering 124 if(center==-3)brvltt=3 ! Face centering 125 if(center==1)brvltt=5 ! A-Face centering 126 if(center==2)brvltt=6 ! B-Face centering 127 if(center==3)brvltt=4 ! C-Face centering 128 if(iholohedry==5)brvltt=7 ! Rhombohedral 129 130 !Produce the symmetry operations, in the axis of the conventional cell 131 nsymconv=nsym 132 if(center/=0)nsymconv=2*nsymconv 133 if(center==-3)nsymconv=4*nsym 134 ABI_ALLOCATE(symrelconv,(3,3,nsymconv)) 135 ABI_ALLOCATE(tnonsconv,(3,nsymconv)) 136 137 !Produce symrel and tnons in conventional axes, 138 !name them symrelconv and tnonsconv 139 rprimdconv(:,1)=bravais(3:5) 140 rprimdconv(:,2)=bravais(6:8) 141 rprimdconv(:,3)=bravais(9:11) 142 143 if(center/=0)rprimdconv(:,:)=rprimdconv(:,:)*half 144 145 axes(:,:)=zero 146 axes(1,1)=one ; axes(2,2)=one ; axes(3,3)=one 147 symrelconv(:,:,1:nsym)=symrel(:,:,1:nsym) 148 !Note that the number of symmetry operations is still nsym 149 call symrelrot(nsym,rprimdconv,axes,symrelconv,tolsym) 150 call xred2xcart(nsym,rprimdconv,tnonsconv,tnons) 151 !Gives the associated translation, with components in the 152 !interval ]-0.5,0.5] . 153 tnonsconv(:,1:nsym)=tnonsconv(:,1:nsym)-nint(tnonsconv(:,1:nsym)-tol6) 154 155 !If the Bravais lattice is centered, duplicate or quadruplicate 156 !the number of symmetry operations, using the Bravais 157 !lattice shifts 158 nshift=1 159 if(center/=0)nshift=2 160 if(center==-3)nshift=4 161 ABI_ALLOCATE(shift,(3,nshift)) 162 shift(:,1)=zero 163 if(center/=0 .and. center/=-3)then 164 shift(:,2)=half 165 if(center==1)shift(1,2)=zero 166 if(center==2)shift(2,2)=zero 167 if(center==3)shift(3,2)=zero 168 else if(center==-3)then 169 shift(:,2)=half ; shift(1,2)=zero 170 shift(:,3)=half ; shift(2,3)=zero 171 shift(:,4)=half ; shift(3,4)=zero 172 end if ! center/=0 or -3 173 if(nshift/=1)then 174 do ishift=2,nshift 175 symrelconv(:,:,(ishift-1)*nsym+1:ishift*nsym)=symrelconv(:,:,1:nsym) 176 do isym=1,nsym 177 tnonsconv(:,(ishift-1)*nsym+isym)=tnonsconv(:,isym)+shift(:,ishift) 178 end do 179 end do ! ishift 180 end if ! nshift/=1 181 182 !At this stage, all the symmetry operations are available, 183 !expressed in the conventional axis, and also include 184 !the Bravais lattive translations, and associated operations... 185 186 n_axes(:)=0 187 188 ABI_ALLOCATE(determinant,(nsymconv)) 189 190 !Get the determinant 191 call symdet(determinant,nsymconv,symrelconv) 192 193 !Get the order of each the symmetry operation, as well as the maximal order 194 !Also, examine whether each symmetry operation is the inversion, or a root 195 !of the inversion (like -3) 196 !Decide which kind of point symmetry operation it is 197 !Finally assign tnonsconv order and decide the space symmetry operation 198 199 ABI_ALLOCATE(t_axes,(nsymconv)) 200 201 do isym=1,nsymconv 202 203 call symcharac(center, determinant(isym), iholohedry, isym, label, & 204 symrelconv(:,:,isym), tnonsconv(:,isym), t_axes(isym)) 205 if (t_axes(isym) == -1) then 206 write(message, '(a,a,i3,a,3(a,3i4,a),a,3es22.12,a,a,3es22.12)' )ch10,& 207 & ' symspgr : problem with isym=',isym,ch10,& 208 & ' symrelconv(:,1,isym)=',symrelconv(:,1,isym),ch10,& 209 & ' symrelconv(:,2,isym)=',symrelconv(:,2,isym),ch10,& 210 & ' symrelconv(:,3,isym)=',symrelconv(:,3,isym),ch10,& 211 & ' tnonsconv(:,isym)=',tnonsconv(:,isym),ch10,& 212 & ' trialt(:)=',trialt(:) 213 call wrtout(std_out,message,'COLL') 214 write(message, '(a,i4,2a)' )& 215 & 'The space symmetry operation number',isym,ch10,& 216 & 'is not a (translated) root of unity' 217 MSG_BUG(message) 218 else if (t_axes(isym) == -2) then 219 write(message, '(a,i0,a)' )'The symmetry operation number ',isym,' is not a root of unity' 220 MSG_BUG(message) 221 end if 222 223 n_axes(t_axes(isym))=n_axes(t_axes(isym))+1 224 225 end do ! isym=1,nsymconv 226 227 if (sum(n_axes)-nsymconv/=0) then 228 write(message, '(7a)' )& 229 & 'Not all the symmetries have been recognized. ',ch10,& 230 & 'This might be due either to an error in the input file',ch10,& 231 & 'or to a BUG in ABINIT',ch10,& 232 & 'Please contact the ABINIT group.' 233 MSG_WARNING(message) 234 end if 235 236 !DEBUG 237 !write(std_out,*)' symspgr : brvltt,nsymconv=',brvltt,nsymconv 238 !write(std_out,*)' n_axes(1:10)=',n_axes(1:10) 239 !write(std_out,*)' n_axes(11:20)=',n_axes(11:20) 240 !write(std_out,*)' n_axes(21:31)=',n_axes(21:31) 241 !ENDDEBUG 242 243 !Treat cases in which the space group cannot be identified on the 244 !basis of n_axes one need additional informations 245 if(brvltt==1)then 246 ! If the bravais lattice is primitive 247 if(nsymconv==4)then 248 n_axest=(/0,0,0,0,0,0,0,1,1,0, 0,0,0,0,0,2,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0/) 249 if(sum((n_axes-n_axest)**2)==0)then ! Spgroup 27 (Pcc2) or 32 (Pba2) 250 write(std_out,*)' symspgr: 27 or 32' 251 additional_info=2 252 ! Select binary axis 253 do isym=1,nsymconv 254 if(t_axes(isym)==8)then 255 ! Find direction of binary axis 256 if(symrelconv(1,1,isym)==1)direction=1 257 if(symrelconv(2,2,isym)==1)direction=2 258 if(symrelconv(3,3,isym)==1)direction=3 259 end if 260 end do 261 ! Examine the projection of the translation vector of the a, b or c mirror planes 262 ! onto the binary axis 263 do isym=1,nsymconv 264 if(t_axes(isym)==16)then 265 if(abs(tnonsconv(direction,isym))>tol8)additional_info=1 266 end if 267 end do 268 end if 269 else if(nsymconv==8)then 270 n_axest=(/0,0,0,0,1,0,0,1,1,0, 0,0,0,0,1,2,0,0,0,2, 0,0,0,0,0,0,0,0,0,0,0/) 271 if(sum((n_axes-n_axest)**2)==0)then ! Spgroup 55 (Pbam) or 57 (Pbcm) 272 write(std_out,*)' symspgr: 55 or 57' 273 additional_info=1 274 ! Select mirror plane m 275 do isym=1,nsymconv 276 if(t_axes(isym)==15)then 277 ! Find direction of mirror plane 278 if(symrelconv(1,1,isym)==-1)direction=1 279 if(symrelconv(2,2,isym)==-1)direction=2 280 if(symrelconv(3,3,isym)==-1)direction=3 281 end if 282 end do 283 ! Examine the projection of the translation vector of the a, b, or c mirror planes 284 ! onto the binary axis 285 do isym=1,nsymconv 286 if(t_axes(isym)==16)then 287 if(abs(tnonsconv(direction,isym))>tol8)additional_info=2 288 end if 289 end do 290 end if 291 n_axest=(/0,0,0,0,1,0,0,1,1,0, 0,0,0,0,0,2,0,1,0,2, 0,0,0,0,0,0,0,0,0,0,0/) 292 if(sum((n_axes-n_axest)**2)==0)then ! Spgroup 56 (Pccn) or 60 (Pbcn) 293 write(std_out,*)' symspgr: 56 or 60' 294 additional_info=1 295 ! Select mirror plane n 296 do isym=1,nsymconv 297 if(t_axes(isym)==18)then 298 ! Find direction of mirror plane 299 if(symrelconv(1,1,isym)==-1)direction=1 300 if(symrelconv(2,2,isym)==-1)direction=2 301 if(symrelconv(3,3,isym)==-1)direction=3 302 end if 303 end do 304 ! Examine the projection of the translation vector of the a, b, or c mirror planes 305 ! onto the binary axis 306 do isym=1,nsymconv 307 if(t_axes(isym)==16)then 308 if(abs(tnonsconv(direction,isym))<tol8)additional_info=2 309 end if 310 end do 311 end if 312 end if 313 else if(brvltt==2)then 314 ! In the few next lines, use additional_info as a flag 315 additional_info=0 316 ! If the bravais lattice is inner-centered 317 if(nsymconv==8)then 318 ! Test spgroup 23 (I222) or 24 (I2_{1}2_{1}2_{1}) 319 n_axest=(/0,0,0,0,0,0,1,1,3,0, 0,0,0,0,0,0,0,0,0,3, 0,0,0,0,0,0,0,0,0,0,0/) 320 if(sum((n_axes-n_axest)**2)==0) additional_info=1 321 else if(nsymconv==24)then 322 ! Test spgroup 197 (I23) or 199 (I2_{1}3) 323 n_axest=(/0,0,0,0,0,0,1,1,3,16, 0,0,0,0,0,0,0,0,0,3, 0,0,0,0,0,0,0,0,0,0,0/) 324 if(sum((n_axes-n_axest)**2)==0) additional_info=1 325 end if 326 if(additional_info==1)then 327 write(std_out,*)' symspgr: (23 or 24) or (197 or 199)' 328 ! Select the three binary axes (they might be 2 or 2_1 !) 329 test_direction(:)=0 330 do isym=1,nsymconv 331 if(t_axes(isym)==20)then 332 ! Find direction of axis 333 do direction=1,3 334 if(symrelconv(direction,direction,isym)==1)then 335 test_direction(direction)=1 336 if(abs(tnonsconv(direction,isym))<tol8)then 337 vect(:,direction)=tnonsconv(:,isym) 338 else 339 vect(:,direction)=tnonsconv(:,isym)+half 340 end if 341 vect(:,direction)=vect(:,direction)-nint(vect(:,direction)-tol8) 342 vect(direction,direction)=zero 343 end if 344 end do ! direction=1,3 345 end if ! if binary axis 346 end do ! isym 347 if(test_direction(1)/=1 .or. test_direction(2)/=1 .and. test_direction(3)/=1)then 348 write(message, '(5a,3i4)' )& 349 & 'For space groups 23, 24, 197 or 197, the three binary axes',ch10,& 350 & 'are not equally partitioned along the x, y and z directions',ch10,& 351 & 'test_direction(1:3)=',test_direction(:) 352 MSG_BUG(message) 353 end if 354 additional_info=1 355 if(abs(vect(1,2)-vect(1,3))>tol8 .or. & 356 & abs(vect(2,1)-vect(2,3))>tol8 .or. & 357 & abs(vect(3,1)-vect(3,2))>tol8) additional_info=2 358 end if ! additional informations are needed 359 end if ! brvltt==1 360 361 if (brvltt==0 .or. brvltt==1) then ! Primitive 362 call symlist_prim(additional_info,nsymconv,n_axes,spgroup) 363 else if(brvltt==2)then 364 call symlist_bcc(additional_info,nsymconv,n_axes,spgroup) 365 else if(brvltt==3)then 366 call symlist_fcc(nsymconv,n_axes,spgroup) 367 else 368 call symlist_others(brvltt,nsymconv,n_axes,spgroup) 369 end if 370 371 if(spgroup==0) then 372 write(message, '(a,a,a,a,a)' )& 373 & 'Could not find the space group.',ch10,& 374 & 'This often happens when the user selects a restricted set of symmetries ',ch10,& 375 & 'in the input file, instead of letting the code automatically find symmetries.' 376 MSG_WARNING(message) 377 end if 378 379 spgorig=1 ; spgaxor=1 380 call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig) 381 382 if(spgroup/=0)then 383 write(message, '(a,i4,2x,a,a,a,a,a)' ) ' symspgr : spgroup=',spgroup,trim(brvsb),trim(intsb),' (=',trim(schsb),')' 384 call wrtout(std_out,message,'COLL') 385 end if 386 387 if(bravais(1)==7)then 388 write(message, '(a)' ) ' symspgr : optical characteristics = isotropic ' 389 call wrtout(std_out,message,'COLL') 390 else if(bravais(1)==4 .or. bravais(1)==5 .or. bravais(1)==6)then 391 write(message, '(a)' ) ' symspgr : optical characteristics = uniaxial ' 392 call wrtout(std_out,message,'COLL') 393 ! Identify the first symmetry operation that is order 3, 4 or 6 394 found=0 395 do isym=1,nsym 396 ! Proper rotations 397 if( minval( abs( t_axes(isym)-(/10,12,14,22,23,24,25,26,27,28,29,30,31/) ))==0) then 398 found=1 ; exit 399 ! Improper symmetry operations 400 else if( minval( abs( t_axes(isym)-(/1,2,3/) ))==0) then 401 found=-1 ; exit 402 end if 403 end do 404 if(found==-1 .or. found==1)then 405 symrel_uni=symrel(:,:,isym) 406 if(found==-1)symrel_uni=-symrel_uni 407 ! Now, symrel_uni is a rotation of order 3, 4, 6, for which the axis must be identified 408 ! It is actually the only eigenvector with eigenvalue 1. It can be found by cross products 409 ! Subtract the unit matrix. 410 do ii=1,3 411 symrel_uni(ii,ii)=symrel_uni(ii,ii)-1 412 end do 413 found=0 414 do ii=1,3 415 jj=ii+1 ; if(jj==4)jj=1 416 ! Cross product 417 uniaxis=symrel_uni(ii,:).x.symrel_uni(jj,:) 418 if(sum(uniaxis**2)/=0)then 419 found=1 ; exit 420 end if 421 end do 422 if(found==1)then 423 ! Try to reduce the length, by an integer factor (try only primes 2, 3, 5, 7, 11) 424 prime=(/2,3,5,7,11/) 425 ii=1 426 do while (ii<6) 427 uniaxis_try=uniaxis/prime(ii) 428 if(sum(abs(uniaxis_try*prime(ii)-uniaxis))==0)then 429 uniaxis=uniaxis_try 430 else 431 ii=ii+1 432 end if 433 end do 434 write(message, '(a,3i4)' ) ' Optical axis (in reduced coordinates, real space ) :',uniaxis 435 end if 436 end if 437 if(found==0)then 438 write(message, '(a)' ) ' However, the axis has not been found. Sorry for this.' 439 end if 440 call wrtout(std_out,message,'COLL') 441 end if 442 443 ABI_DEALLOCATE(determinant) 444 ABI_DEALLOCATE(shift) 445 ABI_DEALLOCATE(symrelconv) 446 ABI_DEALLOCATE(tnonsconv) 447 ABI_DEALLOCATE(t_axes) 448 449 DBG_EXIT("COLL") 450 451 end subroutine symspgr