TABLE OF CONTENTS
ABINIT/symsghexa [ Functions ]
NAME
symsghexa
FUNCTION
Yields all the TRIGONAL & HEXAGONAL symmetry operations starting from the space group symbol. according to the International Tables of Crystallography, 1983.
COPYRIGHT
Copyright (C) 1999-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
msym = default number of symmetries nsym = the number of symmetry operations shubnikov= magnetic type of the space group to be generated spgaxor = ossible orientation of the axes system spgorig = the origin choice (1 or 2) for the axes system spgroup = the numeric symbol of the space groups spgroupma= number of the magnetic space group
OUTPUT
brvltt = bravais lattice type, here, only for rhombohedral groups with hexagonal axes symafm(nsym)=(anti)ferromagnetic part of symmetry operations symrel(3,3,nsym) = 3D matrix containg symmetry operations tnons(3,nsym) = 2D matrix containing translations associated
PARENTS
gensymspgr
CHILDREN
bldgrp,spgdata
SOURCE
41 #if defined HAVE_CONFIG_H 42 #include "config.h" 43 #endif 44 45 #include "abi_common.h" 46 47 48 subroutine symsghexa(brvltt,msym,nsym,shubnikov,spgaxor,spgorig,spgroup,& 49 & spgroupma,symafm,symrel,tnons) 50 51 use defs_basis 52 use m_profiling_abi 53 54 !This section has been created automatically by the script Abilint (TD). 55 !Do not modify the following lines by hand. 56 #undef ABI_FUNC 57 #define ABI_FUNC 'symsghexa' 58 use interfaces_41_geometry, except_this_one => symsghexa 59 !End of the abilint section 60 61 implicit none 62 63 !Arguments ------------------------------------ 64 !scalars 65 integer,intent(in) :: msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma 66 integer,intent(inout) :: brvltt !vz_i 67 !arrays 68 integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i 69 real(dp),intent(inout) :: tnons(3,msym) !vz_i 70 71 !Local variables ----------------------------- 72 !scalars 73 integer :: ii,nogen,sporder 74 real(dp),parameter :: fivesixth=5.0d0/6.0d0,twothird=2.0d0/3.0d0 75 character(len=1) :: brvsb 76 character(len=15) :: intsb,ptintsb,ptschsb,schsb 77 character(len=35) :: intsbl 78 !arrays 79 integer :: genm(3,3),genmmp(3,3),genswm(3,3),genswmmm(3,3),genswmmp(3,3) 80 integer :: genswp(3,3) 81 82 !************************************************************************* 83 84 !DEBUG 85 !write(std_out,*) 'symsghexa',spgroup,shubnikov,spgroupma 86 !ENDDEBUG 87 88 !The identity operation belongs to all space groups 89 symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 90 91 !Predefine some generators 92 genswm(:,:)=0 ; genswm(2,1)=1 ; genswm(1,2)=1 ; genswm(3,3)=-1 !reshape((/0,1,0,1,0,0,0,0,-1/),(/3,3/),(/0,0/),(/2,1/) ) 93 genswmmm(:,:)=0 ; genswmmm(2,1)=-1 ; genswmmm(1,2)=-1 ; genswmmm(3,3)=-1 94 !reshape((/0,-1,0,-1,0,0,0,0,-1/),(/3,3/),(/0,0/),(/2,1/) ) 95 genswmmp(:,:)=0 ; genswmmp(2,1)=-1 ; genswmmp(1,2)=-1 ; genswmmp(3,3)=1 96 !reshape((/0,-1,0,-1,0,0,0,0,1/),(/3,3/),(/0,0/),(/2,1/) ) 97 genswp(:,:)=0 ; genswp(2,1)=1 ; genswp(1,2)=1 ; genswp(3,3)=1 98 !reshape((/0,1,0,1,0,0,0,0,1/),(/3,3/),(/0,0/),(/2,1/) ) 99 genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)=1 100 !reshape((/-1,0,0,0,-1,0,0,0,1/),(/3,3/),(/0,0/),(/2,1/) ) 101 102 !Initialize the associated translations matrix to 0 103 do ii=1,nsym 104 tnons(:,ii)= 0.0d0 105 end do 106 nogen=0 107 108 !Default non-magnetic behaviour 109 symafm(1:nsym)=1 110 111 !************************************************************************* 112 113 !Treat TRIGONAL case 114 if(143<=spgroup .and. spgroup<=167)then 115 116 ! The hexagonal axis choice (orientation) is first treated 117 if (spgaxor == 1) then 118 119 ! This matrix is common to ALL trigonal spatial groups in this orientation 120 ! (Note : this is the 3- symmetry operation) 121 symrel(:,:,2)=0 ; symrel(1,1,2)=-1 ; symrel(1,2,2)=1 ; symrel(2,1,2)=-1 ; symrel(3,3,2)=1 122 ! reshape((/-1,1,0,-1,0,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) ) 123 124 ! Assigns the generators to each space group 125 select case (spgroup) 126 case (143,146,147,148) !P3, R3, PB3, RB3 127 symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(2,2,3)=-1 ; symrel(3,3,3)=1 128 ! reshape((/0,-1,0,1,-1,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) ) 129 nogen=0 130 case (144) !P31 131 tnons(:,2)=(/0.d0,0.d0,twothird/) 132 symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(2,2,3)=-1 ; symrel(3,3,3)=1 133 ! reshape((/0,-1,0,1,-1,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) ) 134 tnons(:,3)=(/0.d0,0.d0,third/) 135 nogen=0 136 case (145) !P32 137 tnons(:,2)=(/0.d0,0.d0,third/) 138 symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(2,2,3)=-1 ; symrel(3,3,3)=1 139 ! reshape((/0,-1,0,1,-1,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) ) 140 tnons(:,3)=(/0.d0,0.d0,twothird/) 141 nogen=0 142 case (149) !P312 143 symrel(:,:,3) = genswmmm(:,:) 144 nogen=3 145 case (150,155) !P321, R32 146 symrel(:,:,3) = genswm(:,:) 147 nogen=3 148 case (151) !P3112 149 tnons(:,2)=(/0.d0,0.d0,twothird/) 150 symrel(:,:,3) = genswmmm(:,:) 151 nogen=3 152 case (152) !P3121 153 tnons(:,2)=(/0.d0,0.d0,twothird/) 154 symrel(:,:,3) = genswm(:,:) 155 nogen=3 156 case (153) !P3212 157 tnons(:,2)=(/0.d0,0.d0,third/) 158 symrel(:,:,3) = genswmmm(:,:) 159 nogen=3 160 case (154) !P3221 161 tnons(:,2)=(/0.d0,0.d0,third/) 162 symrel(:,:,3) = genswm(:,:) 163 nogen=3 164 case (156,160,164,166) !P3m1, R3m, PB3m1, RB3m 165 symrel(:,:,3) = genswmmp(:,:) 166 nogen=3 167 case (157,162) !P31m, PB31m 168 symrel(:,:,3) = genswp(:,:) 169 nogen=3 170 case (158,161,165,167) !P3c1, R3c, PB3c1, RB3c 171 symrel(:,:,3) = genswmmp(:,:) 172 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 173 nogen=3 174 case (159,163) !P31c, PB31c 175 symrel(:,:,3) = genswp(:,:) 176 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 177 nogen=3 178 end select 179 180 select case (spgroup) 181 case (146,148,155,160,166,167) 182 brvltt=7 183 end select 184 185 ! Quite simple, because the generator of even order is always the third one. 186 if(shubnikov==3)then 187 select case(spgroupma) 188 case (23,27,31,35,39,43,47,51,55,59,63,67,71,76,77,82,83,88,89,94,95,& 189 & 100,101,106,107) 190 symafm(3)=-1 191 end select 192 end if 193 194 else if (spgaxor == 2) then 195 ! The rhombohedral axis choice (orientation) is now treated 196 write(std_out,*)'rhombohedral axes' 197 ! Assignment of common three-fold rotation 198 symrel(:,:,2)=0 ; symrel(1,3,2)=1 ; symrel(3,2,2)=1 ; symrel(2,1,2)=1 199 ! reshape((/0,0,1,1,0,0,0,1,0/),(/3,3/),(/0,0/),(/2,1/) ) 200 symrel(:,:,3)=0 ; symrel(3,1,3)=1 ; symrel(2,3,3)=1 ; symrel(1,2,3)=1 201 ! reshape((/0,1,0,0,0,1,1,0,0/), (/3,3/), (/0,0/), (/2,1/) ) 202 203 select case (spgroup) 204 case (146,148) !R3 205 case (155,166) !R32, RB3m 206 symrel(:,:,4) = genswmmm(:,:) 207 nogen=4 208 case (160) !R3m 209 symrel(:,:,4) = genswp(:,:) 210 nogen=4 211 case (161,167) !R3c, RB3c 212 symrel(:,:,4) = genswp(:,:) 213 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 214 nogen=4 215 end select 216 217 if(shubnikov==3)then 218 select case(spgroupma) 219 case (47,67,71,99,101,106,107) 220 symafm(4)=-1 221 end select 222 end if 223 224 ! End selection of axis orientation 225 end if 226 227 ! End trigonal groups 228 end if 229 230 !************************************************************************* 231 232 !Treat HEXAGONAL case 233 if(168<=spgroup .and. spgroup<=194)then 234 235 ! This matrix (6) is common to most hexagonal spatial groups, except 174,187,188,189,190 236 symrel(:,:,2)=0 ; symrel(1,1,2)=1 ; symrel(1,2,2)=-1 ; symrel(2,1,2)=1 ; symrel(3,3,2)=1 237 ! reshape((/1,-1,0,1,0,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) ) 238 ! This one (6 bar) is present in the other cases 239 genm(:,:)=0 ; genm(1,2)=-1 ; genm(2,1)=1 ; genm(2,2)=-1 ; genm(3,3)=-1 240 ! reshape((/0,-1,0,1,-1,0,0,0,-1/), (/3,3/), (/0,0/), (/2,1/) ) 241 select case(spgroup) 242 case (168,175) !P6 243 nogen=2 244 case (169) !P61 245 tnons(:,2)=(/0.d0,0.d0,sixth/) 246 nogen=2 247 case (170) !P65 248 tnons(:,2)=(/0.d0,0.d0,fivesixth/) 249 nogen=2 250 case (171) !P62 251 tnons(:,2)=(/0.d0,0.d0,third/) 252 nogen=2 253 case (172) !P64 254 tnons(:,2)=(/0.d0,0.d0,twothird/) 255 nogen=2 256 case (173,176) !P63, P63/m 257 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 258 nogen=2 259 case (174) !PB6 260 symrel(:,:,2) = genm(:,:) 261 nogen=2 262 case (177) !P622 263 symrel(:,:,3) = genswm(:,:) 264 nogen=3 265 case (178) !P6122 266 tnons(:,2)=(/0.d0,0.d0,sixth/) 267 symrel(:,:,3) = genswm(:,:) 268 nogen=3 269 case (179) !P6522 270 tnons(:,2)=(/0.d0,0.d0,fivesixth/) 271 symrel(:,:,3) = genswm(:,:) 272 nogen=3 273 case (180) !P6222 274 tnons(:,2)=(/0.d0,0.d0,third/) 275 symrel(:,:,3) = genswm(:,:) 276 nogen=3 277 case (181) !P6422 278 tnons(:,2)=(/0.d0,0.d0,twothird/) 279 symrel(:,:,3) = genswm(:,:) 280 nogen=3 281 case (182) !P6322 282 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 283 symrel(:,:,3) = genswm(:,:) 284 nogen=3 285 case (183,191) !P6mm, P6/mmm 286 symrel(:,:,3) = genswp(:,:) 287 nogen=3 288 case (184,192) !P6cc, P6/mcc 289 symrel(:,:,3) = genswp(:,:) 290 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 291 nogen=3 292 case (185,193) !P63cm, P63/mcm 293 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 294 symrel(:,:,3) = genswp(:,:) 295 nogen=3 296 case (186,194) !P63mc, P63/mmc 297 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 298 symrel(:,:,3) = genswp(:,:) 299 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 300 nogen=3 301 case (187) !PB6m2 302 symrel(:,:,2)=0 ; symrel(1,2,2)=-1 ; symrel(2,1,2)=1 ; symrel(2,2,2)=-1 ; symrel(3,3,2)=-1 303 ! reshape((/0,-1,0,1,-1,0,0,0,-1/), (/3,3/), (/0,0/), (/2,1/) ) 304 symrel(:,:,3)=0 ; symrel(1,1,3)=-1 ; symrel(1,2,3)=1 ; symrel(2,2,3)=1 ; symrel(3,3,3)=1 305 ! reshape((/-1,1,0,0,1,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) ) 306 nogen=3 307 if (shubnikov==3) then 308 if (spgroupma==211) symafm(2:3)=-1 309 if (spgroupma==212) symafm(2)=-1 310 if (spgroupma==213) symafm(3)=-1 311 end if 312 case (188) !PB6c2 313 symrel(:,:,2) = genm(:,:) 314 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 315 symrel(:,:,3) = genswmmm(:,:) 316 nogen=3 317 case (189) !PB62m 318 symrel(:,:,2) = genm(:,:) 319 symrel(:,:,3) = genswp(:,:) 320 nogen=3 321 case (190) !PB62c 322 symrel(:,:,2) = genm(:,:) 323 symrel(:,:,3) = genswp(:,:) 324 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 325 nogen=3 326 end select 327 328 if(shubnikov==3)then 329 select case(spgroupma) 330 ! spgroup from 168 to 176 are OK, 177 to 194 are not done 331 case (111,115,119,123,127,131,135,139,141,145,147,152,158,164,170,& 332 & 176,182,187,193,199,205,217,224,230,237,239,247,249,257,259,267,269) 333 symafm(2)=-1 334 case(153,159,165,171,177,183,189,195,& 335 & 201,207,219,225,231,240,241,250,251,260,261,270,271) 336 symafm(3)=-1 337 case(151,157,163,169,175,181,188,194,200,206,218,223,229,236,238,246,248,256,258,266,268) 338 symafm(2:3)=-1 339 end select 340 end if 341 342 call spgdata(brvsb,intsb,intsbl,ptintsb,& 343 & ptschsb,schsb,spgaxor,spgroup,sporder,spgorig) 344 345 ! End HEXAGONAL groups 346 end if 347 348 !*************************************************************************** 349 350 !DEBUG 351 !write(std_out,*) 'symsghexa : out with nogen = ',nogen 352 !ENDDEBUG 353 354 355 if (nogen>0) then 356 call bldgrp(msym,nogen,nsym,symafm,symrel,tnons) 357 end if 358 359 !DEBUG 360 !write(std_out,*)'symrel:' 361 !write(std_out,*) symrel(:,:,1:nsym) 362 !ENDDEBUG 363 364 end subroutine symsghexa