TABLE OF CONTENTS
ABINIT/symsgortho [ Functions ]
NAME
symsgortho
FUNCTION
Yields all the ORTHORHOMBIC symmetry operations starting from the space group symbol. It deals only with the orthorhombic groups taken in the standard orientation 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 spgorig = the origin choice (1 or 2) for the axes system spgaxor = the possible orientation of the axes system spgroup = the numeric symbol of the space groups spgroupma= number of the magnetic space group
OUTPUT
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 symsgortho(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 'symsgortho' 58 use interfaces_41_geometry, except_this_one => symsgortho 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 66 integer,intent(in) :: spgroupma 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 ! nogen = number of generators selected 73 !scalars 74 integer :: nogen,sporder 75 character(len=1) :: brvsb 76 character(len=15) :: intsb,ptintsb,ptschsb,schsb 77 character(len=35) :: intsbl 78 !arrays 79 integer :: genmmm(3,3),genmmp(3,3),genmpm(3,3),genmpp(3,3),genpmm(3,3) 80 integer :: genpmp(3,3),genppm(3,3) 81 82 ! ************************************************************************* 83 84 !DEBUG 85 !write(std_out,*)'symsgortho ( orthorhombic groups) : enter with space group ',spgroup 86 !ENDDEBUG 87 88 !The orientation of the space group: 89 !first we will permute the input coordinates of the atoms, xred 90 !then we will make the calculation in the "normal" space group 91 !then the coordinates are reoriented to match the initial orientation 92 !and finally the symrel is reoriented to correspond to the new orientation 93 !further all the calculations are performed into the space group 94 !with the user-defined orientation 95 96 !The identity operation belongs to all space groups 97 symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 98 99 nogen=4 100 101 !Predefine some generators 102 genmpp(:,:)=0 ; genmpp(1,1)=-1 ; genmpp(2,2)= 1 ; genmpp(3,3)= 1 103 genpmp(:,:)=0 ; genpmp(1,1)= 1 ; genpmp(2,2)=-1 ; genpmp(3,3)= 1 104 genppm(:,:)=0 ; genppm(1,1)= 1 ; genppm(2,2)= 1 ; genppm(3,3)=-1 105 genpmm(:,:)=0 ; genpmm(1,1)= 1 ; genpmm(2,2)=-1 ; genpmm(3,3)=-1 106 genmpm(:,:)=0 ; genmpm(1,1)=-1 ; genmpm(2,2)= 1 ; genmpm(3,3)=-1 107 genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)= 1 108 genmmm(:,:)=0 ; genmmm(1,1)=-1 ; genmmm(2,2)=-1 ; genmmm(3,3)=-1 109 110 111 !For all the groups in this routine symrel(:,:,2) is the same 112 symrel(:,:,2)=genmmp(:,:) 113 114 !Default non-magnetic behaviour 115 symafm(1:nsym)=1 116 117 !DEBUG 118 !write(std_out,*) 'symsgortho:',spgroup,shubnikov,spgroupma 119 !ENDDEBUG 120 121 !assigns the generators to each space group 122 select case (spgroup) 123 ! ORTHORHOMBIC space groups 124 case (16,21,22,23) !P222, C222, F222, I222 125 symrel(:,:,3) = genmpm(:,:) 126 symrel(:,:,4) = genpmm(:,:) 127 if(shubnikov==3)symafm(3:4)=-1 128 case (17,20) !P2221, C2221 129 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 130 symrel(:,:,3) = genpmm(:,:) 131 symrel(:,:,4) = genmpm(:,:) 132 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 133 if(shubnikov==3)then 134 symafm(4)=-1 135 if(spgroupma==9) symafm(3)=-1 136 if(spgroupma==10)symafm(2)=-1 137 if(spgroupma==33)symafm(3:4)=-1 138 if(spgroupma==34)symafm(2)=-1 139 end if 140 case (18) !P21212 141 symrel(:,:,3) = genmpm(:,:) 142 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 143 symrel(:,:,4) = genpmm(:,:) 144 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 145 if(shubnikov==3)then 146 symafm(3)=-1 147 if(spgroupma==18)symafm(4)=-1 148 if(spgroupma==19)symafm(2)=-1 149 end if 150 case (19,24) !P212121, I212121 151 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 152 symrel(:,:,3) = genmpm(:,:) 153 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 154 symrel(:,:,4) = genpmm(:,:) 155 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 156 if(shubnikov==3)symafm(3:4)=-1 157 case (25,35,38,42,44) !Pmm2, Cmm2, Amm2, Fmm2, Imm2 158 symrel(:,:,3) = genpmp(:,:) 159 symrel(:,:,4) = genmpp(:,:) 160 case (26,36) !Pmc21, Cmc21 161 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 162 symrel(:,:,3) = genpmp(:,:) 163 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 164 symrel(:,:,4) = genmpp(:,:) 165 case (27,37) !Pcc2, Ccc2 166 symrel(:,:,3) = genpmp(:,:) 167 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 168 symrel(:,:,4) = genmpp(:,:) 169 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 170 case (28,40,46) !Pma2, Ama2, Ima2 171 symrel(:,:,3) = genpmp(:,:) 172 tnons(:,3)=(/0.5d0,0.d0,0.d0/) 173 symrel(:,:,4) = genmpp(:,:) 174 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 175 case (29) !Pca21 176 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 177 symrel(:,:,3) = genpmp(:,:) 178 tnons(:,3)=(/0.5d0,0.d0,0.d0/) 179 symrel(:,:,4) = genmpp(:,:) 180 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 181 case (30) !Pnc2 182 symrel(:,:,3) = genpmp(:,:) 183 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 184 symrel(:,:,4) = genmpp(:,:) 185 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 186 case (31) !Pmn21 187 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 188 symrel(:,:,3) = genpmp(:,:) 189 tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 190 symrel(:,:,4) = genmpp(:,:) 191 case (32,41,45) !Pba2, Aba2, Iba2 192 symrel(:,:,3) = genpmp(:,:) 193 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 194 symrel(:,:,4) = genmpp(:,:) 195 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 196 case (33) !Pna21 197 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 198 symrel(:,:,3) = genpmp(:,:) 199 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 200 symrel(:,:,4) = genmpp(:,:) 201 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 202 case (34) !Pnn2 203 symrel(:,:,3) = genpmp(:,:) 204 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 205 symrel(:,:,4) = genmpp(:,:) 206 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 207 case (39) !Abm2 208 symrel(:,:,3) = genpmp(:,:) 209 tnons(:,3)=(/0.d0,0.5d0,0.d0/) 210 symrel(:,:,4) = genmpp(:,:) 211 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 212 case (43) !Fdd2 213 symrel(:,:,3) = genpmp(:,:) 214 tnons(:,3)=(/0.25d0,0.25d0,0.25d0/) 215 symrel(:,:,4) = genmpp(:,:) 216 tnons(:,4)=(/0.25d0,0.25d0,0.25d0/) 217 case (47,65,69,71) !Pmmm, Cmmm, Fmmm, Immm 218 symrel(:,:,3) = genmpm(:,:) 219 symrel(:,:,4) = genpmm(:,:) 220 case (48) !Pnnn 221 symrel(:,:,3) = genmpm(:,:) 222 symrel(:,:,4) = genpmm(:,:) 223 symrel(:,:,5) = genmmm(:,:) 224 symrel(:,:,6) = genppm(:,:) 225 symrel(:,:,7) = genpmp(:,:) 226 symrel(:,:,8) = genmpp(:,:) 227 if (spgorig==1) then 228 tnons(:,5)=(/0.5d0,0.5d0,0.5d0/) 229 tnons(:,6)=(/0.5d0,0.5d0,0.5d0/) 230 tnons(:,7)=(/0.5d0,0.5d0,0.5d0/) 231 tnons(:,8)=(/0.5d0,0.5d0,0.5d0/) 232 else 233 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 234 tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 235 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 236 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 237 tnons(:,7)=(/0.5d0,0.d0,0.5d0/) 238 tnons(:,8)=(/0.d0,0.5d0,0.5d0/) 239 end if 240 if(shubnikov==3)then 241 if(spgroupma==259)then 242 symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1 243 else if(spgroupma==260)then 244 symafm(3:4)=-1 ; symafm(7:8)=-1 245 else if(spgroupma==261)then 246 symafm(5:8)=-1 247 end if 248 end if 249 nogen=0 250 case (49,66) !Pccm, Cccm 251 symrel(:,:,4) = genpmm(:,:) 252 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 253 symrel(:,:,3) = genmpm(:,:) 254 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 255 case (50) !Pban 256 symrel(:,:,3) = genmpm(:,:) 257 symrel(:,:,4) = genpmm(:,:) 258 symrel(:,:,5) = genmmm(:,:) 259 symrel(:,:,6) = genppm(:,:) 260 symrel(:,:,7) = genpmp(:,:) 261 symrel(:,:,8) = genmpp(:,:) 262 if (spgorig==1) then 263 tnons(:,5)=(/0.5d0,0.5d0,0.d0/) 264 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 265 tnons(:,7)=(/0.5d0,0.5d0,0.d0/) 266 tnons(:,8)=(/0.5d0,0.5d0,0.d0/) 267 else 268 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 269 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 270 tnons(:,3)=(/0.5d0,0.d0,0.d0/) 271 tnons(:,7)=(/0.5d0,0.d0,0.d0/) 272 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 273 tnons(:,8)=(/0.d0,0.5d0,0.d0/) 274 end if 275 if(shubnikov==3)then 276 if(spgroupma==279)then 277 symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1 278 else if(spgroupma==280)then 279 symafm(3:4)=-1 ; symafm(5:6)=-1 280 else if(spgroupma==281)then 281 symafm(3:4)=-1 ; symafm(7:8)=-1 282 else if(spgroupma==282)then 283 symafm(2:8:2)=-1 284 else if(spgroupma==283)then 285 symafm(5:8)=-1 286 end if 287 end if 288 nogen=0 289 case (51) !Pmma 290 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 291 symrel(:,:,3) = genmpm(:,:) 292 symrel(:,:,4) = genpmm(:,:) 293 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 294 case (52) !Pnna 295 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 296 symrel(:,:,3) = genmpm(:,:) 297 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 298 symrel(:,:,4) = genpmm(:,:) 299 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 300 case (53) !Pmna 301 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 302 symrel(:,:,3) = genmpm(:,:) 303 tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 304 symrel(:,:,4) = genpmm(:,:) 305 case (54) !Pcca 306 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 307 symrel(:,:,3) = genmpm(:,:) 308 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 309 symrel(:,:,4) = genpmm(:,:) 310 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 311 case (55,72) !Pbam, Ibam 312 symrel(:,:,3) = genmpm(:,:) 313 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 314 symrel(:,:,4) = genpmm(:,:) 315 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 316 case (56) !Pccn 317 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 318 symrel(:,:,3) = genmpm(:,:) 319 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 320 symrel(:,:,4) = genpmm(:,:) 321 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 322 case (57) !Pbcm 323 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 324 symrel(:,:,3) = genmpm(:,:) 325 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 326 symrel(:,:,4) = genpmm(:,:) 327 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 328 case (58) !Pnnm 329 symrel(:,:,3) = genmpm(:,:) 330 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 331 symrel(:,:,4) = genpmm(:,:) 332 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 333 case (59) !Pmmn 334 symrel(:,:,3) = genmpm(:,:) 335 symrel(:,:,4) = genpmm(:,:) 336 symrel(:,:,5) = genmmm(:,:) 337 symrel(:,:,6) = genppm(:,:) 338 symrel(:,:,7) = genpmp(:,:) 339 symrel(:,:,8) = genmpp(:,:) 340 if (spgorig==1) then 341 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 342 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 343 tnons(:,5)=(/0.5d0,0.5d0,0.d0/) 344 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 345 else 346 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 347 tnons(:,3)=(/0.d0,0.5d0,0.d0/) 348 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 349 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 350 tnons(:,7)=(/0.d0,0.5d0,0.d0/) 351 tnons(:,8)=(/0.5d0,0.d0,0.d0/) 352 end if 353 if(shubnikov==3)then 354 if(spgroupma==407)then 355 symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1 356 else if(spgroupma==408)then 357 symafm(3:4)=-1 ; symafm(5:6)=-1 358 else if(spgroupma==409)then 359 symafm(3:4)=-1 ; symafm(7:8)=-1 360 else if(spgroupma==410)then 361 symafm(2:8:2)=-1 362 else if(spgroupma==411)then 363 symafm(5:8)=-1 364 end if 365 end if 366 nogen=0 367 case (60) !Pbcn 368 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 369 symrel(:,:,3) = genmpm(:,:) 370 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 371 symrel(:,:,4) = genpmm(:,:) 372 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 373 case (61,73) !Pbca, Ibca 374 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 375 symrel(:,:,3) = genmpm(:,:) 376 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 377 symrel(:,:,4) = genpmm(:,:) 378 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 379 case (62) !Pnma 380 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 381 symrel(:,:,3) = genmpm(:,:) 382 tnons(:,3)=(/0.d0,0.5d0,0.d0/) 383 symrel(:,:,4) = genpmm(:,:) 384 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 385 case (63) !Cmcm 386 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 387 symrel(:,:,3) = genmpm(:,:) 388 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 389 symrel(:,:,4) = genpmm(:,:) 390 if(shubnikov==3)then 391 if(spgroupma==459 .or. spgroupma==463)symafm(2:3)=-1 392 if(spgroupma==460 .or. spgroupma==464)symafm(4)=-1 393 if(spgroupma==460 .or. spgroupma==464)symafm(2)=-1 394 if(spgroupma==461 .or. spgroupma==462)symafm(3:4)=-1 395 end if 396 case (64) !Cmca 397 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 398 symrel(:,:,3) = genmpm(:,:) 399 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 400 symrel(:,:,4) = genpmm(:,:) 401 case (67) !Cmma 402 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 403 symrel(:,:,4) = genpmm(:,:) 404 symrel(:,:,3) = genmpm(:,:) 405 tnons(:,3)=(/0.d0,0.5d0,0.d0/) 406 case (68) !Ccca 407 symrel(:,:,3) = genmpm(:,:) 408 symrel(:,:,4) = genpmm(:,:) 409 symrel(:,:,5) = genmmm(:,:) 410 symrel(:,:,6) = genppm(:,:) 411 symrel(:,:,7) = genpmp(:,:) 412 symrel(:,:,8) = genmpp(:,:) 413 if (spgorig==1) then 414 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 415 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 416 tnons(:,5)=(/0.d0,0.5d0,0.5d0/) 417 tnons(:,6)=(/0.5d0,0.d0,0.5d0/) 418 tnons(:,7)=(/0.d0,0.5d0,0.5d0/) 419 tnons(:,8)=(/0.5d0,0.d0,0.5d0/) 420 else 421 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 422 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 423 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 424 tnons(:,6)=(/0.5d0,0.d0,0.d0/) 425 tnons(:,7)=(/0.d0,0.d0,0.5d0/) 426 tnons(:,8)=(/0.5d0,0.d0,0.5d0/) 427 end if 428 if(shubnikov==3)then 429 if(spgroupma==513)then 430 symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1 431 else if(spgroupma==514)then 432 symafm(3:4)=-1 ; symafm(5:6)=-1 433 else if(spgroupma==515)then 434 symafm(3:4)=-1 ; symafm(7:8)=-1 435 else if(spgroupma==516)then 436 symafm(2:8:2)=-1 437 else if(spgroupma==517)then 438 symafm(5:8)=-1 439 end if 440 end if 441 nogen=0 442 case (70) !Fddd 443 symrel(:,:,3) = genmpm(:,:) 444 symrel(:,:,4) = genpmm(:,:) 445 symrel(:,:,5) = genmmm(:,:) 446 symrel(:,:,6) = genppm(:,:) 447 symrel(:,:,7) = genpmp(:,:) 448 symrel(:,:,8) = genmpp(:,:) 449 if (spgorig==1) then 450 tnons(:,5)=(/0.25d0,0.25d0,0.25d0/) 451 tnons(:,6)=(/0.25d0,0.25d0,0.25d0/) 452 tnons(:,7)=(/0.25d0,0.25d0,0.25d0/) 453 tnons(:,8)=(/0.25d0,0.25d0,0.25d0/) 454 else 455 tnons(:,2)=(/0.75d0,0.75d0,0.d0/) 456 tnons(:,3)=(/0.75d0,0.d0,0.75d0/) 457 tnons(:,4)=(/0.d0,0.75d0,0.75d0/) 458 ! JWZ DEBUGGING BEGIN 459 ! the following lines were present in 5.7 as of Dec 10 2008 but 460 ! gave wrong results in a spgroup 70 spgorig 2 case (Na2SO4) 461 ! tnons(:,5)=(/0.25d0,0.25d0,0.d0/) ! original code 462 ! tnons(:,6)=(/0.25d0,0.d0,0.25d0/) ! original code 463 ! tnons(:,7)=(/0.d0,0.25d0,0.25d0/) ! original code 464 ! here are the corrected values of tnons for this case 465 tnons(:,5)=(/0.0d0,0.0d0,0.0d0/) 466 tnons(:,6)=(/0.25d0,0.25d0,0.0d0/) 467 tnons(:,7)=(/0.25d0,0.0d0,0.25d0/) 468 tnons(:,8)=(/0.0d0,0.25d0,0.25d0/) 469 ! JWZ DEBUGGING END 470 end if 471 if(shubnikov==3)then 472 if(spgroupma==529)then 473 symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1 474 else if(spgroupma==530)then 475 symafm(3:4)=-1 ; symafm(7:8)=-1 476 else if(spgroupma==531)then 477 symafm(5:8)=-1 478 end if 479 end if 480 nogen=0 481 case (74) !Imma 482 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 483 symrel(:,:,3) = genmpm(:,:) 484 tnons(:,3)=(/0.d0,0.5d0,0.d0/) 485 symrel(:,:,4) = genpmm(:,:) 486 end select 487 488 if (shubnikov==3) then 489 select case (spgroupma) 490 case (59,68,80,89,101,113,125,137,146,158,167,174,182,189,197,& 491 & 205,213,221,226,231,237,243,270,292,296,308,312,324,328,340,344,& 492 & 358,370,380,384,420,424,444,448,460,464,472,476) 493 symafm(2)=-1 494 symafm(4)=-1 495 case (69,90,102,114,126,147,175,190,198,206,214,244) 496 symafm(2:3)=-1 497 case (60,70,81,91,103,115,127,138,148,159,168,176,183,191,199,& 498 & 207,215,222,227,232,238,245,252,268,269,293,294,309,310,& 499 & 325,326,341,342,356,357,368,369,381,382,396,397,421,422,436,445,& 500 & 446,461,462,473,474,484,485,494,495,504,505,524,536,& 501 & 542,543,551,557,558) 502 symafm(3:4)=-1 503 case (251,267,291,295,307,311,323,327,339,343,355,367,379,383,& 504 & 395,398,419,423,435,443,447,459,463,& 505 & 471,475,483,486,493,496,503,506,523,535,541,544,550,556,559) 506 symafm(2:3)=-1 507 508 end select 509 end if 510 511 if (nogen>1)then 512 call bldgrp(msym,nogen,nsym,symafm,symrel,tnons) 513 end if 514 515 call spgdata(brvsb,intsb,intsbl,ptintsb,& 516 & ptschsb,schsb,spgaxor,spgroup,sporder,spgorig) 517 518 !DEBUG 519 !write(std_out,*)'symsgortho : end of symmetry assignement' 520 !ENDDEBUG 521 522 end subroutine symsgortho