TABLE OF CONTENTS
ABINIT/symsgmono [ Functions ]
NAME
symsgmono
FUNCTION
Yields all the MONOCLINIC symmetry operations starting from the space group symbol. according to the International Tables of Crystallography, 1983. It solves also the problem of the axes orientation according to the spgaxor
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 = the orientation choice of the unit cell spgorig = possible origin of 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 (1=P; 2=I; 3=F; 4=C; 5=A; 6=B; 7=R) 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
spgdata
SOURCE
43 #if defined HAVE_CONFIG_H 44 #include "config.h" 45 #endif 46 47 #include "abi_common.h" 48 49 50 subroutine symsgmono(brvltt,msym,nsym,shubnikov,spgaxor,spgorig,spgroup,& 51 & spgroupma,symafm,symrel,tnons) 52 53 54 use defs_basis 55 use m_profiling_abi 56 57 !This section has been created automatically by the script Abilint (TD). 58 !Do not modify the following lines by hand. 59 #undef ABI_FUNC 60 #define ABI_FUNC 'symsgmono' 61 use interfaces_41_geometry, except_this_one => symsgmono 62 !End of the abilint section 63 64 implicit none 65 66 !Arguments ------------------------------------ 67 !scalars 68 integer,intent(in) :: msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma 69 integer,intent(inout) :: brvltt !vz_i 70 !arrays 71 integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i 72 real(dp),intent(inout) :: tnons(3,msym) !vz_i 73 74 !Local variables ----------------------------- 75 !scalars 76 integer :: sporder 77 character(len=1) :: brvsb 78 character(len=15) :: intsb,ptintsb,ptschsb,schsb 79 character(len=35) :: intsbl 80 !arrays 81 integer :: genmmm(3,3),genmmp(3,3),genmpm(3,3),genmpp(3,3),genpmm(3,3) 82 integer :: genpmp(3,3),genppm(3,3) 83 84 ! ************************************************************************* 85 !the identity operation belonging to all space groups 86 symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 87 88 !Predefine some generators 89 genmpp(:,:)=0 ; genmpp(1,1)=-1 ; genmpp(2,2)= 1 ; genmpp(3,3)= 1 90 genpmp(:,:)=0 ; genpmp(1,1)= 1 ; genpmp(2,2)=-1 ; genpmp(3,3)= 1 91 genppm(:,:)=0 ; genppm(1,1)= 1 ; genppm(2,2)= 1 ; genppm(3,3)=-1 92 genpmm(:,:)=0 ; genpmm(1,1)= 1 ; genpmm(2,2)=-1 ; genpmm(3,3)=-1 93 genmpm(:,:)=0 ; genmpm(1,1)=-1 ; genmpm(2,2)= 1 ; genmpm(3,3)=-1 94 genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)= 1 95 genmmm(:,:)=0 ; genmmm(1,1)=-1 ; genmmm(2,2)=-1 ; genmmm(3,3)=-1 96 97 !Default non-magnetic behaviour 98 symafm(1:nsym)=1 99 100 !assigns the generators to each space group 101 select case (spgroup) 102 case (3) ! P2 103 select case (spgaxor) 104 case (1) ! 3:b, P2_b = P2 105 symrel(:,:,2) = genmpm(:,:) 106 case (2) ! 3:a, P2_a = P2 107 symrel(:,:,2) = genpmm(:,:) 108 case (3) ! 3:c, P2_c = P2 109 symrel(:,:,2) = genmmp(:,:) 110 end select 111 if(shubnikov==3)symafm(2)=-1 112 case (4) ! P21 113 select case (spgaxor) 114 case (1) ! 3:b, P21_b = P21 115 symrel(:,:,2) = genmpm(:,:) 116 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 117 case (2) ! 3:a, P21_a = P21 118 symrel(:,:,2) = genpmm(:,:) 119 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 120 case (3) ! 3:c, P21_c = P21 121 symrel(:,:,2) = genmmp(:,:) 122 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 123 end select 124 if(shubnikov==3)symafm(2)=-1 125 case (5) ! C2 126 select case (spgaxor) 127 case (1) ! 5:b1, C2 = C2 128 symrel(:,:,2) = genmpm(:,:) 129 brvltt=4 130 case (2) ! 5:a1, B2_a = C2 131 symrel(:,:,2) = genpmm(:,:) 132 brvltt=6 133 case (3) ! 5:a2, C2_a = C2 134 symrel(:,:,2) = genpmm(:,:) 135 brvltt=4 136 case (4) ! 5:a3, I2_a = C2 137 symrel(:,:,2) = genpmm(:,:) 138 brvltt=2 139 case (5) ! 5:b2, A2_b = C2 140 symrel(:,:,2) = genmpm(:,:) 141 brvltt=5 142 case (6) ! 5:b3, I2_b = C2 143 symrel(:,:,2) = genmpm(:,:) 144 brvltt=2 145 case (7) ! 5:c1, A2_c = C2 146 symrel(:,:,2) = genmmp(:,:) 147 brvltt=5 148 case (8) ! 5:c2, B2_c = C2 149 symrel(:,:,2) = genmmp(:,:) 150 brvltt=6 151 case (9) ! 5:c3, I2_c = C2 152 symrel(:,:,2) = genmmp(:,:) 153 brvltt=2 154 end select 155 if(shubnikov==3)symafm(2)=-1 156 case (6) ! Pm 157 select case (spgaxor) 158 case (1) ! Pm_b = Pm 159 symrel(:,:,2) = genpmp(:,:) 160 case (2) ! Pm_a = Pm 161 symrel(:,:,2) = genmpp(:,:) 162 case (3) ! Pm_c = Pm 163 symrel(:,:,2) = genppm(:,:) 164 end select 165 if(shubnikov==3)symafm(2)=-1 166 case (7) ! Pc 167 select case (spgaxor) 168 case (1) ! 7:b1, Pc_b = Pc 169 symrel(:,:,2) = genpmp(:,:) 170 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 171 case (2) ! 7:a1, Pb_a = Pc 172 symrel(:,:,2) = genmpp(:,:) 173 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 174 case (3) ! 7:a2, Pn_a = Pc 175 symrel(:,:,2) = genmpp(:,:) 176 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 177 case (4) ! 7:a3, Pc_a = Pc 178 symrel(:,:,2) = genmpp(:,:) 179 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 180 case (5) ! 7:b2, Pn_b = Pc 181 symrel(:,:,2) = genpmp(:,:) 182 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 183 case (6) ! 7:b3, Pa_b = Pc 184 symrel(:,:,2) = genpmp(:,:) 185 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 186 case (7) ! 7:c1, Pa_c = Pc 187 symrel(:,:,2) = genppm(:,:) 188 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 189 case (8) ! 7:c2, Pn_c = Pc 190 symrel(:,:,2) = genppm(:,:) 191 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 192 case (9) ! 7:c3, Pb_c = Pb 193 symrel(:,:,2) = genppm(:,:) 194 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 195 end select 196 if(shubnikov==3)symafm(2)=-1 197 case (8) ! Cm 198 select case (spgaxor) 199 case (1) ! 8:b1, Cm = Cm 200 symrel(:,:,2) = genpmp(:,:) 201 brvltt=4 202 case (2) ! 8:a1, Bm_a = Cm 203 symrel(:,:,2) = genmpp(:,:) 204 brvltt=6 205 case (3) ! 8:a2, Cm_a = Cm 206 symrel(:,:,2) = genmpp(:,:) 207 brvltt=4 208 case (4) ! 8:a3, Im_a = Cm 209 symrel(:,:,2) = genmpp(:,:) 210 brvltt=2 211 case (5) ! 8:b2, Am_b = Cm 212 symrel(:,:,2) = genpmp(:,:) 213 brvltt=5 214 case (6) ! 8:b3, Im_b = Cm 215 symrel(:,:,2) = genpmp(:,:) 216 brvltt=2 217 case (7) ! 8:c1, Am_c = Cm 218 symrel(:,:,2) = genppm(:,:) 219 brvltt=5 220 case (8) ! 8:c2, Bm_c = Bm 221 symrel(:,:,2) = genppm(:,:) 222 brvltt=6 223 case (9) ! 8:c3, Im_c = Cm 224 symrel(:,:,2) = genppm(:,:) 225 brvltt=2 226 end select 227 if(shubnikov==3)symafm(2)=-1 228 case (9) ! Cc 229 select case (spgaxor) 230 case (1) ! 9:b1, Cc_b = Cc 231 symrel(:,:,2) = genpmp(:,:) 232 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 233 brvltt=4 234 case (2) ! 9:a1, Bb_a = Cc 235 symrel(:,:,2) = genmpp(:,:) 236 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 237 brvltt=6 238 case (3) ! 9:a2, Cn_a = Cc 239 symrel(:,:,2) = genmpp(:,:) 240 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 241 brvltt=4 242 case (4) ! 9:a3, Ic_a = Cc 243 symrel(:,:,2) = genmpp(:,:) 244 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 245 brvltt=2 246 case (5) ! 9:b2, An_b = Cc 247 symrel(:,:,2) = genpmp(:,:) 248 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 249 brvltt=5 250 case (6) ! 9:b3, Ia_b = Cc 251 symrel(:,:,2) = genpmp(:,:) 252 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 253 brvltt=2 254 case (7) ! 9:c1, Aa_c = Cc 255 symrel(:,:,2) = genppm(:,:) 256 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 257 brvltt=5 258 case (8) ! 9:c2, B(b+c)_c = Cc 259 symrel(:,:,2) = genppm(:,:) 260 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 261 brvltt=6 262 case (9) ! 9:c3, Ib_c = Cc 263 symrel(:,:,2) = genppm(:,:) 264 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 265 brvltt=2 266 end select 267 if(shubnikov==3)symafm(2)=-1 268 case (10) ! P2/m 269 select case (spgaxor) 270 case (1) ! 10:b, P2/m = P2/m 271 symrel(:,:,2) = genmpm(:,:) 272 symrel(:,:,3) = genmmm(:,:) 273 symrel(:,:,4) = genpmp(:,:) 274 case (2) ! 10:a, P2/m_a = P2/m 275 symrel(:,:,2) = genpmm(:,:) 276 symrel(:,:,3) = genmmm(:,:) 277 symrel(:,:,4) = genmpp(:,:) 278 case (3) ! 10:c, P2/m_c = P2/m 279 symrel(:,:,2) = genmmp(:,:) 280 symrel(:,:,3) = genmmm(:,:) 281 symrel(:,:,4) = genppm(:,:) 282 end select 283 if(shubnikov==3)then 284 symafm(2:4)=-1 ! Default 285 if(spgroupma==44)symafm(4)=1 286 if(spgroupma==45)symafm(2)=1 287 if(spgroupma==46)symafm(3)=1 288 end if 289 case (11) ! P21/m 290 select case (spgaxor) 291 case (1) ! 11:b, P21/m = P21/m 292 symrel(:,:,2) = genmpm(:,:) 293 symrel(:,:,3) = genmmm(:,:) 294 symrel(:,:,4) = genpmp(:,:) 295 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 296 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 297 case (2) ! 11:a, P21/m_a = P21/m 298 symrel(:,:,2) = genpmm(:,:) 299 symrel(:,:,3) = genmmm(:,:) 300 symrel(:,:,4) = genmpp(:,:) 301 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 302 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 303 case (3) ! 11:c, P21/m_c = P21/m 304 symrel(:,:,2) = genmmp(:,:) 305 symrel(:,:,3) = genmmm(:,:) 306 symrel(:,:,4) = genppm(:,:) 307 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 308 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 309 end select 310 if(shubnikov==3)then 311 symafm(2:4)=-1 ! Default 312 if(spgroupma==52)symafm(4)=1 313 if(spgroupma==53)symafm(2)=1 314 if(spgroupma==54)symafm(3)=1 315 end if 316 case (12) ! C2/m 317 select case (spgaxor) 318 case (1) ! 12:b1, C2/m = C2/m 319 symrel(:,:,2) = genmpm(:,:) 320 symrel(:,:,3) = genmmm(:,:) 321 symrel(:,:,4) = genpmp(:,:) 322 brvltt=4 323 case (2) ! 12:a1, B2/m_a = C2/m 324 symrel(:,:,2) = genpmm(:,:) 325 symrel(:,:,3) = genmmm(:,:) 326 symrel(:,:,4) = genmpp(:,:) 327 brvltt=6 328 case (3) ! 12:a2, C2/m_a = C2/m 329 symrel(:,:,2) = genpmm(:,:) 330 symrel(:,:,3) = genmmm(:,:) 331 symrel(:,:,4) = genmpp(:,:) 332 brvltt=4 333 case (4) ! 12:a3, I2/m_a = C2/m 334 symrel(:,:,2) = genpmm(:,:) 335 symrel(:,:,3) = genmmm(:,:) 336 symrel(:,:,4) = genmpp(:,:) 337 brvltt=2 338 case (5) ! 12:b2, A2/m_b = C2/m 339 symrel(:,:,2) = genmpm(:,:) 340 symrel(:,:,3) = genmmm(:,:) 341 symrel(:,:,4) = genpmp(:,:) 342 brvltt=5 343 case (6) ! 12:b3, I2/m_b = C2/m 344 symrel(:,:,2) = genmpm(:,:) 345 symrel(:,:,3) = genmmm(:,:) 346 symrel(:,:,4) = genpmp(:,:) 347 brvltt=2 348 case (7) ! 12:c1, A2/m_c = C2/m 349 symrel(:,:,2) = genmmp(:,:) 350 symrel(:,:,3) = genmmm(:,:) 351 symrel(:,:,4) = genppm(:,:) 352 brvltt=5 353 case (8) ! 12:c2, B2/m_c = B2/m 354 symrel(:,:,2) = genmmp(:,:) 355 symrel(:,:,3) = genmmm(:,:) 356 symrel(:,:,4) = genppm(:,:) 357 brvltt=6 358 case (9) ! 12:c3, I2/m_c = C2/m 359 symrel(:,:,2) = genmmp(:,:) 360 symrel(:,:,3) = genmmm(:,:) 361 symrel(:,:,4) = genppm(:,:) 362 brvltt=2 363 end select 364 if(shubnikov==3)then 365 symafm(2:4)=-1 ! Default 366 if(spgroupma==60)symafm(4)=1 367 if(spgroupma==61)symafm(2)=1 368 if(spgroupma==62)symafm(3)=1 369 end if 370 case (13) ! P2/c 371 select case (spgaxor) 372 case (1) ! 13:b1, P2/c = P2/c 373 symrel(:,:,2) = genmpm(:,:) 374 symrel(:,:,3) = genmmm(:,:) 375 symrel(:,:,4) = genpmp(:,:) 376 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 377 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 378 case (2) ! 13:a1, P2/b_a = P2/c 379 symrel(:,:,2) = genpmm(:,:) 380 symrel(:,:,3) = genmmm(:,:) 381 symrel(:,:,4) = genmpp(:,:) 382 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 383 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 384 case (3) ! 13:a2, P2/n_a = P2/c 385 symrel(:,:,2) = genpmm(:,:) 386 symrel(:,:,3) = genmmm(:,:) 387 symrel(:,:,4) = genmpp(:,:) 388 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 389 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 390 case (4) ! 13:a3, P2/c_a = P2/c 391 symrel(:,:,2) = genpmm(:,:) 392 symrel(:,:,3) = genmmm(:,:) 393 symrel(:,:,4) = genmpp(:,:) 394 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 395 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 396 case (5) ! 13:b2, P2/n_b = P2/c 397 symrel(:,:,2) = genmpm(:,:) 398 symrel(:,:,3) = genmmm(:,:) 399 symrel(:,:,4) = genpmp(:,:) 400 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 401 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 402 case (6) ! 13:b3, P2/a_b = P2/c 403 symrel(:,:,2) = genmpm(:,:) 404 symrel(:,:,3) = genmmm(:,:) 405 symrel(:,:,4) = genpmp(:,:) 406 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 407 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 408 case (7) ! 13:c1, P2/a_c = P2/c 409 symrel(:,:,2) = genmmp(:,:) 410 symrel(:,:,3) = genmmm(:,:) 411 symrel(:,:,4) = genppm(:,:) 412 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 413 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 414 case (8) ! 13:c2, P2/n_c = P2/c 415 symrel(:,:,2) = genmmp(:,:) 416 symrel(:,:,3) = genmmm(:,:) 417 symrel(:,:,4) = genppm(:,:) 418 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 419 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 420 case (9) ! 13:c3, P2/b_c = P2/b 421 symrel(:,:,2) = genmmp(:,:) 422 symrel(:,:,3) = genmmm(:,:) 423 symrel(:,:,4) = genppm(:,:) 424 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 425 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 426 end select 427 if(shubnikov==3)then 428 symafm(2:4)=-1 ! Default 429 if(spgroupma==67)symafm(4)=1 430 if(spgroupma==68)symafm(2)=1 431 if(spgroupma==69)symafm(3)=1 432 end if 433 case (14) ! P21/c 434 select case (spgaxor) 435 case (1) ! 14:b1, P21/c_b = P21/c 436 symrel(:,:,2) = genmpm(:,:) 437 symrel(:,:,3) = genmmm(:,:) 438 symrel(:,:,4) = genpmp(:,:) 439 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 440 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 441 case (2) ! 14:a1, P21/a_b = P21/c 442 symrel(:,:,2) = genpmm(:,:) 443 symrel(:,:,3) = genmmm(:,:) 444 symrel(:,:,4) = genmpp(:,:) 445 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 446 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 447 case (3) ! 14:a2, P21/n_a = P21/c 448 symrel(:,:,2) = genpmm(:,:) 449 symrel(:,:,3) = genmmm(:,:) 450 symrel(:,:,4) = genmpp(:,:) 451 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 452 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 453 case (4) ! 14:a3, P21/c_a = P21/c 454 symrel(:,:,2) = genpmm(:,:) 455 symrel(:,:,3) = genmmm(:,:) 456 symrel(:,:,4) = genmpp(:,:) 457 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 458 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 459 case (5) ! 14:b2, P21/n_b = P21/c 460 symrel(:,:,2) = genmpm(:,:) 461 symrel(:,:,3) = genmmm(:,:) 462 symrel(:,:,4) = genpmp(:,:) 463 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 464 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 465 case (6) ! 14:b3, P21/a_b = P21/c 466 symrel(:,:,2) = genmpm(:,:) 467 symrel(:,:,3) = genmmm(:,:) 468 symrel(:,:,4) = genpmp(:,:) 469 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 470 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 471 case (7) ! 14:c1, P21/a_c = P21/c 472 symrel(:,:,2) = genmmp(:,:) 473 symrel(:,:,3) = genmmm(:,:) 474 symrel(:,:,4) = genppm(:,:) 475 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 476 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 477 case (8) ! 14:c2, P21/n_c = P21/c 478 symrel(:,:,2) = genmmp(:,:) 479 symrel(:,:,3) = genmmm(:,:) 480 symrel(:,:,4) = genppm(:,:) 481 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 482 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 483 case (9) ! 14/c3, P21/b_c = P21/b 484 symrel(:,:,2) = genmmp(:,:) 485 symrel(:,:,3) = genmmm(:,:) 486 symrel(:,:,4) = genppm(:,:) 487 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 488 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 489 end select 490 if(shubnikov==3)then 491 symafm(2:4)=-1 ! Default 492 if(spgroupma==77)symafm(4)=1 493 if(spgroupma==78)symafm(2)=1 494 if(spgroupma==79)symafm(3)=1 495 end if 496 case (15) ! C2/c 497 select case (spgaxor) 498 case (1) ! 15:b1, C2/c_b = C2/c 499 symrel(:,:,2) = genmpm(:,:) 500 symrel(:,:,3) = genmmm(:,:) 501 symrel(:,:,4) = genpmp(:,:) 502 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 503 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 504 brvltt = 4 505 case (2) ! 15:a1, B2/b_a = C2/c 506 symrel(:,:,2) = genpmm(:,:) 507 symrel(:,:,3) = genmmm(:,:) 508 symrel(:,:,4) = genmpp(:,:) 509 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 510 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 511 brvltt = 6 512 case (3) ! 15:a2, C2/n_a = C2/c 513 symrel(:,:,2) = genpmm(:,:) 514 symrel(:,:,3) = genmmm(:,:) 515 symrel(:,:,4) = genmpp(:,:) 516 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 517 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 518 brvltt = 4 519 case (4) ! 15:a3, I2/c_a = C2/c 520 symrel(:,:,2) = genpmm(:,:) 521 symrel(:,:,3) = genmmm(:,:) 522 symrel(:,:,4) = genmpp(:,:) 523 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 524 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 525 brvltt = 2 526 case (5) ! 15:b2, A2/n_b = C2/c 527 symrel(:,:,2) = genmpm(:,:) 528 symrel(:,:,3) = genmmm(:,:) 529 symrel(:,:,4) = genpmp(:,:) 530 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 531 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 532 brvltt = 5 533 case (6) ! 15:b3, I2/a_b = C2/c 534 symrel(:,:,2) = genmpm(:,:) 535 symrel(:,:,3) = genmmm(:,:) 536 symrel(:,:,4) = genpmp(:,:) 537 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 538 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 539 brvltt = 2 540 case (7) ! 15:c1, A2/a_c = C2/c 541 symrel(:,:,2) = genmmp(:,:) 542 symrel(:,:,3) = genmmm(:,:) 543 symrel(:,:,4) = genppm(:,:) 544 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 545 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 546 brvltt = 5 547 case (8) ! 15:c2, B21/b_c = C2/c 548 symrel(:,:,2) = genmmp(:,:) 549 symrel(:,:,3) = genmmm(:,:) 550 symrel(:,:,4) = genppm(:,:) 551 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 552 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 553 brvltt = 6 554 case (9) ! 15:c3, I2/b_c = C2/c 555 symrel(:,:,2) = genmmp(:,:) 556 symrel(:,:,3) = genmmm(:,:) 557 symrel(:,:,4) = genppm(:,:) 558 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 559 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 560 brvltt = 2 561 end select 562 if(shubnikov==3)then 563 symafm(2:4)=-1 ! Default 564 if(spgroupma==87)symafm(4)=1 565 if(spgroupma==88)symafm(2)=1 566 if(spgroupma==89)symafm(3)=1 567 end if 568 end select 569 570 call spgdata(brvsb,intsb,intsbl,ptintsb,& 571 & ptschsb,schsb,spgaxor,spgroup,sporder,spgorig) 572 573 end subroutine symsgmono