TABLE OF CONTENTS
ABINIT/gensymshub [ Functions ]
NAME
gensymshub
FUNCTION
Analyse the Shubnikov space group, from the input of the the Fedorov space group number, and the Shubnikov space group number: 1) determine the type (III or IV) 2) for type (IV), determine the translation generating the anti-ferromagnetic operations At present, follow strictly the specifications of Bradley and Cracknell. However, one should take into account the orientation of the symmetry group (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
spgroup = number of space group spgroupma = number of magnetic space group
OUTPUT
genafm(3) = in case of shubnikov type IV, translation, generator of the anti-ferromagnetic symmetry operations shubnikov = type of the shubnikov group
PARENTS
ingeo
CHILDREN
SOURCE
39 #if defined HAVE_CONFIG_H 40 #include "config.h" 41 #endif 42 43 #include "abi_common.h" 44 45 46 subroutine gensymshub(genafm,spgroup,spgroupma,shubnikov) 47 48 use defs_basis 49 use m_errors 50 use m_profiling_abi 51 52 !This section has been created automatically by the script Abilint (TD). 53 !Do not modify the following lines by hand. 54 #undef ABI_FUNC 55 #define ABI_FUNC 'gensymshub' 56 !End of the abilint section 57 58 implicit none 59 60 !Arguments ------------------------------------ 61 !scalars 62 integer,intent(in) :: spgroup,spgroupma 63 integer,intent(out) :: shubnikov 64 !arrays 65 real(dp),intent(out) :: genafm(3) 66 67 !Local variables ------------------------------ 68 !scalars 69 integer :: brvlttbw=0,spgrmatch=1 70 character(len=500) :: message 71 72 ! ************************************************************************* 73 74 !List of the input parameters 75 !DEBUG 76 !write(std_out,*) 77 !write(std_out,*)' gensymshub : enter with:' 78 !write(std_out,*)' brvlttbw = ',brvlttbw 79 !write(std_out,*)' spgroup = ',spgroup 80 !write(std_out,*)' spgroupma = ',spgroupma 81 !ENDDEBUG 82 83 !Test for consistency the magnetic and non-magnetic space group 84 select case (spgroup) 85 case(1) 86 if (.not.(spgroupma>=3 .and. 3 >=spgroupma) ) spgrmatch=0 87 case(2) 88 if (.not.(spgroupma>=6 .and. 7 >=spgroupma) ) spgrmatch=0 89 case(3) 90 if (.not.(spgroupma>=3 .and. 6 >=spgroupma) ) spgrmatch=0 91 case(4) 92 if (.not.(spgroupma>=9 .and. 12 >=spgroupma) ) spgrmatch=0 93 case(5) 94 if (.not.(spgroupma>=15 .and. 17 >=spgroupma) ) spgrmatch=0 95 case(6) 96 if (.not.(spgroupma>=20 .and. 23 >=spgroupma) ) spgrmatch=0 97 case(7) 98 if (.not.(spgroupma>=26 .and. 31 >=spgroupma) ) spgrmatch=0 99 case(8) 100 if (.not.(spgroupma>=34 .and. 36 >=spgroupma) ) spgrmatch=0 101 case(9) 102 if (.not.(spgroupma>=39 .and. 41 >=spgroupma) ) spgrmatch=0 103 case(10) 104 if (.not.(spgroupma>=44 .and. 49 >=spgroupma) ) spgrmatch=0 105 case(11) 106 if (.not.(spgroupma>=52 .and. 57 >=spgroupma) ) spgrmatch=0 107 case(12) 108 if (.not.(spgroupma>=60 .and. 64 >=spgroupma) ) spgrmatch=0 109 case(13) 110 if (.not.(spgroupma>=67 .and. 74 >=spgroupma) ) spgrmatch=0 111 case(14) 112 if (.not.(spgroupma>=77 .and. 84 >=spgroupma) ) spgrmatch=0 113 case(15) 114 if (.not.(spgroupma>=87 .and. 91 >=spgroupma) ) spgrmatch=0 115 case(16) 116 if (.not.(spgroupma>= 3 .and. 6>=spgroupma) ) spgrmatch=0 117 case(17) 118 if (.not.(spgroupma>= 9 .and. 15 >=spgroupma) ) spgrmatch=0 119 case(18) 120 if (.not.(spgroupma>=18 .and. 24 >=spgroupma) ) spgrmatch=0 121 case(19) 122 if (.not.(spgroupma>=27 .and. 30 >=spgroupma) ) spgrmatch=0 123 case(20) 124 if (.not.(spgroupma>=33 .and. 37 >=spgroupma) ) spgrmatch=0 125 case(21) 126 if (.not.(spgroupma>=40 .and. 44 >=spgroupma) ) spgrmatch=0 127 case(22) 128 if (.not.(spgroupma>=47 .and. 48 >=spgroupma) ) spgrmatch=0 129 case(23) 130 if (.not.(spgroupma>=51 .and. 52 >=spgroupma) ) spgrmatch=0 131 case(24) 132 if (.not.(spgroupma>=55 .and. 56 >=spgroupma) ) spgrmatch=0 133 case(25) 134 if (.not.(spgroupma>=59 .and. 65 >=spgroupma) ) spgrmatch=0 135 case(26) 136 if (.not.(spgroupma>=68 .and. 77 >=spgroupma) ) spgrmatch=0 137 case(27) 138 if (.not.(spgroupma>=80 .and. 86 >=spgroupma) ) spgrmatch=0 139 case(28) 140 if (.not.(spgroupma>=89 .and. 98 >=spgroupma) ) spgrmatch=0 141 case(29) 142 if (.not.(spgroupma>=101 .and. 110 >=spgroupma) ) spgrmatch=0 143 case(30) 144 if (.not.(spgroupma>=113 .and. 122 >=spgroupma) ) spgrmatch=0 145 case(31) 146 if (.not.(spgroupma>=125 .and. 134 >=spgroupma) ) spgrmatch=0 147 case(32) 148 if (.not.(spgroupma>=137 .and. 143 >=spgroupma) ) spgrmatch=0 149 case(33) 150 if (.not.(spgroupma>=146 .and. 155 >=spgroupma) ) spgrmatch=0 151 case(34) 152 if (.not.(spgroupma>=158 .and. 164 >=spgroupma) ) spgrmatch=0 153 case(35) 154 if (.not.(spgroupma>=167 .and. 171>=spgroupma) ) spgrmatch=0 155 case(36) 156 if (.not.(spgroupma>=174 .and. 179>=spgroupma) ) spgrmatch=0 157 case(37) 158 if (.not.(spgroupma>=182 .and. 186 >=spgroupma) ) spgrmatch=0 159 case(38) 160 if (.not.(spgroupma>=189 .and. 194 >=spgroupma) ) spgrmatch=0 161 case(39) 162 if (.not.(spgroupma>=197 .and. 202 >=spgroupma) ) spgrmatch=0 163 case(40) 164 if (.not.(spgroupma>=205 .and. 210 >=spgroupma) ) spgrmatch=0 165 case(41) 166 if (.not.(spgroupma>=213 .and. 218 >=spgroupma) ) spgrmatch=0 167 case(42) 168 if (.not.(spgroupma>=221 .and. 223>=spgroupma) ) spgrmatch=0 169 case(43) 170 if (.not.(spgroupma>=226 .and. 228>=spgroupma) ) spgrmatch=0 171 case(44) 172 if (.not.(spgroupma>=231 .and. 234 >=spgroupma) ) spgrmatch=0 173 case(45) 174 if (.not.(spgroupma>=237 .and. 240 >=spgroupma) ) spgrmatch=0 175 case(46) 176 if (.not.(spgroupma>=243 .and. 248 >=spgroupma) ) spgrmatch=0 177 case(47) 178 if (.not.(spgroupma>=251 .and. 256 >=spgroupma) ) spgrmatch=0 179 case(48) 180 if (.not.(spgroupma>=259 .and. 264 >=spgroupma) ) spgrmatch=0 181 case(49) 182 if (.not.(spgroupma>=267 .and. 276 >=spgroupma) ) spgrmatch=0 183 case(50) 184 if (.not.(spgroupma>=279 .and. 288 >=spgroupma) ) spgrmatch=0 185 case(51) 186 if (.not.(spgroupma>=291 .and. 304 >=spgroupma) ) spgrmatch=0 187 case(52) 188 if (.not.(spgroupma>=307 .and. 320 >=spgroupma) ) spgrmatch=0 189 case(53) 190 if (.not.(spgroupma>=323 .and. 336 >=spgroupma) ) spgrmatch=0 191 case(54) 192 if (.not.(spgroupma>=339 .and. 352 >=spgroupma) ) spgrmatch=0 193 case(55) 194 if (.not.(spgroupma>=355 .and. 364 >=spgroupma) ) spgrmatch=0 195 case(56) 196 if (.not.(spgroupma>=367 .and. 376 >=spgroupma) ) spgrmatch=0 197 case(57) 198 if (.not.(spgroupma>=379 .and. 392 >=spgroupma) ) spgrmatch=0 199 case(58) 200 if (.not.(spgroupma>=395 .and. 404 >=spgroupma) ) spgrmatch=0 201 case(59) 202 if (.not.(spgroupma>=407 .and. 416 >=spgroupma) ) spgrmatch=0 203 case(60) 204 if (.not.(spgroupma>=419 .and. 432 >=spgroupma) ) spgrmatch=0 205 case(61) 206 if (.not.(spgroupma>=435 .and. 440 >=spgroupma) ) spgrmatch=0 207 case(62) 208 if (.not.(spgroupma>=443 .and. 456 >=spgroupma) ) spgrmatch=0 209 case(63) 210 if (.not.(spgroupma>=459 .and. 468 >=spgroupma) ) spgrmatch=0 211 case(64) 212 if (.not.(spgroupma>=471 .and. 480 >=spgroupma) ) spgrmatch=0 213 case(65) 214 if (.not.(spgroupma>=483 .and. 490 >=spgroupma) ) spgrmatch=0 215 case(66) 216 if (.not.(spgroupma>=493 .and. 500 >=spgroupma) ) spgrmatch=0 217 case(67) 218 if (.not.(spgroupma>=503 .and. 510 >=spgroupma) ) spgrmatch=0 219 case(68) 220 if (.not.(spgroupma>=513 .and. 520 >=spgroupma) ) spgrmatch=0 221 case(69) 222 if (.not.(spgroupma>=523 .and. 526 >=spgroupma) ) spgrmatch=0 223 case(70) 224 if (.not.(spgroupma>=529 .and. 532 >=spgroupma) ) spgrmatch=0 225 case(71) 226 if (.not.(spgroupma>=535 .and. 538 >=spgroupma) ) spgrmatch=0 227 case(72) 228 if (.not.(spgroupma>=541 .and. 547 >=spgroupma) ) spgrmatch=0 229 case(73) 230 if (.not.(spgroupma>=550 .and. 553 >=spgroupma) ) spgrmatch=0 231 case(74) 232 if (.not.(spgroupma>=556 .and. 562 >=spgroupma) ) spgrmatch=0 233 case(75) 234 if (.not.(spgroupma>= 3 .and. 6>=spgroupma) ) spgrmatch=0 235 case(76) 236 if (.not.(spgroupma>= 9 .and. 12 >=spgroupma) ) spgrmatch=0 237 case(77) 238 if (.not.(spgroupma>= 15 .and. 18>=spgroupma) ) spgrmatch=0 239 case(78) 240 if (.not.(spgroupma>= 21 .and. 24>=spgroupma) ) spgrmatch=0 241 case(79) 242 if (.not.(spgroupma>= 27 .and. 28 >=spgroupma) ) spgrmatch=0 243 case(80) 244 if (.not.(spgroupma>= 31 .and. 32 >=spgroupma) ) spgrmatch=0 245 case(81) 246 if (.not.(spgroupma>= 35 .and. 38 >=spgroupma) ) spgrmatch=0 247 case(82) 248 if (.not.(spgroupma>= 41 .and. 42 >=spgroupma) ) spgrmatch=0 249 case(83) 250 if (.not.(spgroupma>=45 .and. 50>=spgroupma) ) spgrmatch=0 251 case(84) 252 if (.not.(spgroupma>=53 .and. 58>=spgroupma) ) spgrmatch=0 253 case(85) 254 if (.not.(spgroupma>=61 .and. 66>=spgroupma) ) spgrmatch=0 255 case(86) 256 if (.not.(spgroupma>=69 .and. 74>=spgroupma) ) spgrmatch=0 257 case(87) 258 if (.not.(spgroupma>=77 .and. 80>=spgroupma) ) spgrmatch=0 259 case(88) 260 if (.not.(spgroupma>=83 .and. 86>=spgroupma) ) spgrmatch=0 261 case(89) 262 if (.not.(spgroupma>=89 .and. 94>=spgroupma) ) spgrmatch=0 263 case(90) 264 if (.not.(spgroupma>=97 .and. 102>=spgroupma) ) spgrmatch=0 265 case(91) 266 if (.not.(spgroupma>=105 .and. 110>=spgroupma) ) spgrmatch=0 267 case(92) 268 if (.not.(spgroupma>=113 .and. 118>=spgroupma) ) spgrmatch=0 269 case(93) 270 if (.not.(spgroupma>=121 .and. 126>=spgroupma) ) spgrmatch=0 271 case(94) 272 if (.not.(spgroupma>=129 .and. 134>=spgroupma) ) spgrmatch=0 273 case(95) 274 if (.not.(spgroupma>=137 .and. 142>=spgroupma) ) spgrmatch=0 275 case(96) 276 if (.not.(spgroupma>=145 .and. 150>=spgroupma) ) spgrmatch=0 277 case(97) 278 if (.not.(spgroupma>=153 .and. 156>=spgroupma) ) spgrmatch=0 279 case(98) 280 if (.not.(spgroupma>=159 .and. 162>=spgroupma) ) spgrmatch=0 281 case(99) 282 if (.not.(spgroupma>=165 .and. 170>=spgroupma) ) spgrmatch=0 283 case(100) 284 if (.not.(spgroupma>=173 .and. 178>=spgroupma) ) spgrmatch=0 285 case(101) 286 if (.not.(spgroupma>=181 .and. 186>=spgroupma) ) spgrmatch=0 287 case(102) 288 if (.not.(spgroupma>=189 .and. 194>=spgroupma) ) spgrmatch=0 289 case(103) 290 if (.not.(spgroupma>=197 .and. 202>=spgroupma) ) spgrmatch=0 291 case(104) 292 if (.not.(spgroupma>=205 .and. 210>=spgroupma) ) spgrmatch=0 293 case(105) 294 if (.not.(spgroupma>=213 .and. 218>=spgroupma) ) spgrmatch=0 295 case(106) 296 if (.not.(spgroupma>=221 .and. 226>=spgroupma) ) spgrmatch=0 297 case(107) 298 if (.not.(spgroupma>=229 .and. 232>=spgroupma) ) spgrmatch=0 299 case(108) 300 if (.not.(spgroupma>=235 .and. 238>=spgroupma) ) spgrmatch=0 301 case(109) 302 if (.not.(spgroupma>=241 .and. 244>=spgroupma) ) spgrmatch=0 303 case(110) 304 if (.not.(spgroupma>=247 .and. 250>=spgroupma) ) spgrmatch=0 305 case(111) 306 if (.not.(spgroupma>=253 .and. 258>=spgroupma) ) spgrmatch=0 307 case(112) 308 if (.not.(spgroupma>=261 .and. 266>=spgroupma) ) spgrmatch=0 309 case(113) 310 if (.not.(spgroupma>=269 .and. 274>=spgroupma) ) spgrmatch=0 311 case(114) 312 if (.not.(spgroupma>=277 .and. 282>=spgroupma) ) spgrmatch=0 313 case(115) 314 if (.not.(spgroupma>=285 .and. 290>=spgroupma) ) spgrmatch=0 315 case(116) 316 if (.not.(spgroupma>=293 .and. 298>=spgroupma) ) spgrmatch=0 317 case(117) 318 if (.not.(spgroupma>=301 .and. 306>=spgroupma) ) spgrmatch=0 319 case(118) 320 if (.not.(spgroupma>=309 .and. 314>=spgroupma) ) spgrmatch=0 321 case(119) 322 if (.not.(spgroupma>=317 .and. 320>=spgroupma) ) spgrmatch=0 323 case(120) 324 if (.not.(spgroupma>=323 .and. 326>=spgroupma) ) spgrmatch=0 325 case(121) 326 if (.not.(spgroupma>=329 .and. 332>=spgroupma) ) spgrmatch=0 327 case(122) 328 if (.not.(spgroupma>=335 .and. 338>=spgroupma) ) spgrmatch=0 329 case(123) 330 if (.not.(spgroupma>=341 .and. 350>=spgroupma) ) spgrmatch=0 331 case(124) 332 if (.not.(spgroupma>=353 .and. 362>=spgroupma) ) spgrmatch=0 333 case(125) 334 if (.not.(spgroupma>=365 .and. 374>=spgroupma) ) spgrmatch=0 335 case(126) 336 if (.not.(spgroupma>=377 .and. 386>=spgroupma) ) spgrmatch=0 337 case(127) 338 if (.not.(spgroupma>=389 .and. 398>=spgroupma) ) spgrmatch=0 339 case(128) 340 if (.not.(spgroupma>=401 .and. 410>=spgroupma) ) spgrmatch=0 341 case(129) 342 if (.not.(spgroupma>=413 .and. 422>=spgroupma) ) spgrmatch=0 343 case(130) 344 if (.not.(spgroupma>=425 .and. 434>=spgroupma) ) spgrmatch=0 345 case(131) 346 if (.not.(spgroupma>=437 .and. 446>=spgroupma) ) spgrmatch=0 347 case(132) 348 if (.not.(spgroupma>=449 .and. 458>=spgroupma) ) spgrmatch=0 349 case(133) 350 if (.not.(spgroupma>=461 .and. 470>=spgroupma) ) spgrmatch=0 351 case(134) 352 if (.not.(spgroupma>=473 .and. 482>=spgroupma) ) spgrmatch=0 353 case(135) 354 if (.not.(spgroupma>=485 .and. 494>=spgroupma) ) spgrmatch=0 355 case(136) 356 if (.not.(spgroupma>=497 .and. 506>=spgroupma) ) spgrmatch=0 357 case(137) 358 if (.not.(spgroupma>=509 .and. 518>=spgroupma) ) spgrmatch=0 359 case(138) 360 if (.not.(spgroupma>=521 .and. 530>=spgroupma) ) spgrmatch=0 361 case(139) 362 if (.not.(spgroupma>=533 .and. 540>=spgroupma) ) spgrmatch=0 363 case(140) 364 if (.not.(spgroupma>=543 .and. 550>=spgroupma) ) spgrmatch=0 365 case(141) 366 if (.not.(spgroupma>=553 .and. 560>=spgroupma) ) spgrmatch=0 367 case(142) 368 if (.not.(spgroupma>=563 .and. 570>=spgroupma) ) spgrmatch=0 369 case(143) 370 if (.not.(spgroupma>=3 .and. 3>=spgroupma) ) spgrmatch=0 371 case(144) 372 if (.not.(spgroupma>=6 .and. 6>=spgroupma) ) spgrmatch=0 373 case(145) 374 if (.not.(spgroupma>=9 .and. 9>=spgroupma) ) spgrmatch=0 375 case(146) 376 if (.not.(spgroupma>=12 .and. 12>=spgroupma) ) spgrmatch=0 377 case(147) 378 if (.not.(spgroupma>=15 .and. 16>=spgroupma) ) spgrmatch=0 379 case(148) 380 if (.not.(spgroupma>=19 .and. 20>=spgroupma) ) spgrmatch=0 381 case(149) 382 if (.not.(spgroupma>=23 .and. 24>=spgroupma) ) spgrmatch=0 383 case(150) 384 if (.not.(spgroupma>=27 .and. 28>=spgroupma) ) spgrmatch=0 385 case(151) 386 if (.not.(spgroupma>=31 .and. 32>=spgroupma) ) spgrmatch=0 387 case(152) 388 if (.not.(spgroupma>=35 .and. 36>=spgroupma) ) spgrmatch=0 389 case(153) 390 if (.not.(spgroupma>=39 .and. 40>=spgroupma) ) spgrmatch=0 391 case(154) 392 if (.not.(spgroupma>=43 .and. 44>=spgroupma) ) spgrmatch=0 393 case(155) 394 if (.not.(spgroupma>=47 .and. 48>=spgroupma) ) spgrmatch=0 395 case(156) 396 if (.not.(spgroupma>=51 .and. 52>=spgroupma) ) spgrmatch=0 397 case(157) 398 if (.not.(spgroupma>=55 .and. 56>=spgroupma) ) spgrmatch=0 399 case(158) 400 if (.not.(spgroupma>=59 .and. 60>=spgroupma) ) spgrmatch=0 401 case(159) 402 if (.not.(spgroupma>=63 .and. 64>=spgroupma) ) spgrmatch=0 403 case(160) 404 if (.not.(spgroupma>=67 .and. 68>=spgroupma) ) spgrmatch=0 405 case(161) 406 if (.not.(spgroupma>=71 .and. 72>=spgroupma) ) spgrmatch=0 407 case(162) 408 if (.not.(spgroupma>=75 .and. 78>=spgroupma) ) spgrmatch=0 409 case(163) 410 if (.not.(spgroupma>=81 .and. 84>=spgroupma) ) spgrmatch=0 411 case(164) 412 if (.not.(spgroupma>=87 .and. 90>=spgroupma) ) spgrmatch=0 413 case(165) 414 if (.not.(spgroupma>=93 .and. 96>=spgroupma) ) spgrmatch=0 415 case(166) 416 if (.not.(spgroupma>=99 .and. 102>=spgroupma) ) spgrmatch=0 417 case(167) 418 if (.not.(spgroupma>=105 .and. 108>=spgroupma) ) spgrmatch=0 419 case(168) 420 if (.not.(spgroupma>=111 .and. 112>=spgroupma) ) spgrmatch=0 421 case(169) 422 if (.not.(spgroupma>=115 .and. 116>=spgroupma) ) spgrmatch=0 423 case(170) 424 if (.not.(spgroupma>=119 .and. 120>=spgroupma) ) spgrmatch=0 425 case(171) 426 if (.not.(spgroupma>=123 .and. 124>=spgroupma) ) spgrmatch=0 427 case(172) 428 if (.not.(spgroupma>=127 .and. 128>=spgroupma) ) spgrmatch=0 429 case(173) 430 if (.not.(spgroupma>=131 .and. 132>=spgroupma) ) spgrmatch=0 431 case(174) 432 if (.not.(spgroupma>=135 .and. 136>=spgroupma) ) spgrmatch=0 433 case(175) 434 if (.not.(spgroupma>=139 .and. 142>=spgroupma) ) spgrmatch=0 435 case(176) 436 if (.not.(spgroupma>=145 .and. 148>=spgroupma) ) spgrmatch=0 437 case(177) 438 if (.not.(spgroupma>=151 .and. 154>=spgroupma) ) spgrmatch=0 439 case(178) 440 if (.not.(spgroupma>=157 .and. 160>=spgroupma) ) spgrmatch=0 441 case(179) 442 if (.not.(spgroupma>=163 .and. 166>=spgroupma) ) spgrmatch=0 443 case(180) 444 if (.not.(spgroupma>=169 .and. 172>=spgroupma) ) spgrmatch=0 445 case(181) 446 if (.not.(spgroupma>=175 .and. 178>=spgroupma) ) spgrmatch=0 447 case(182) 448 if (.not.(spgroupma>=181 .and. 184>=spgroupma) ) spgrmatch=0 449 case(183) 450 if (.not.(spgroupma>=187 .and. 190>=spgroupma) ) spgrmatch=0 451 case(184) 452 if (.not.(spgroupma>=193 .and. 196>=spgroupma) ) spgrmatch=0 453 case(185) 454 if (.not.(spgroupma>=199 .and. 202>=spgroupma) ) spgrmatch=0 455 case(186) 456 if (.not.(spgroupma>=205 .and. 208>=spgroupma) ) spgrmatch=0 457 case(187) 458 if (.not.(spgroupma>=211 .and. 214>=spgroupma) ) spgrmatch=0 459 case(188) 460 if (.not.(spgroupma>=217 .and. 220>=spgroupma) ) spgrmatch=0 461 case(189) 462 if (.not.(spgroupma>=223 .and. 226>=spgroupma) ) spgrmatch=0 463 case(190) 464 if (.not.(spgroupma>=229 .and. 232>=spgroupma) ) spgrmatch=0 465 case(191) 466 if (.not.(spgroupma>=235 .and. 242>=spgroupma) ) spgrmatch=0 467 case(192) 468 if (.not.(spgroupma>=245 .and. 252>=spgroupma) ) spgrmatch=0 469 case(193) 470 if (.not.(spgroupma>=255 .and. 262>=spgroupma) ) spgrmatch=0 471 case(194) 472 if (.not.(spgroupma>=265 .and. 272>=spgroupma) ) spgrmatch=0 473 case(195) 474 if (.not.(spgroupma>=3 .and. 3>=spgroupma) ) spgrmatch=0 475 case(196) 476 if (.not.(spgroupma>=6 .and. 6>=spgroupma) ) spgrmatch=0 477 case(197) 478 spgrmatch=0 479 case(198) 480 if (.not.(spgroupma>=11 .and. 11>=spgroupma) ) spgrmatch=0 481 case(199) 482 spgrmatch=0 483 case(200) 484 if (.not.(spgroupma>=16 .and. 17>=spgroupma) ) spgrmatch=0 485 case(201) 486 if (.not.(spgroupma>=20 .and. 21>=spgroupma) ) spgrmatch=0 487 case(202) 488 if (.not.(spgroupma>=24 .and. 25>=spgroupma) ) spgrmatch=0 489 case(203) 490 if (.not.(spgroupma>=28 .and. 29>=spgroupma) ) spgrmatch=0 491 case(204) 492 if (.not.(spgroupma>=32 .and. 32>=spgroupma) ) spgrmatch=0 493 case(205) 494 if (.not.(spgroupma>=35 .and. 36>=spgroupma) ) spgrmatch=0 495 case(206) 496 if (.not.(spgroupma>=39 .and. 39>=spgroupma) ) spgrmatch=0 497 case(207) 498 if (.not.(spgroupma>=42 .and. 43>=spgroupma) ) spgrmatch=0 499 case(208) 500 if (.not.(spgroupma>=46 .and. 47>=spgroupma) ) spgrmatch=0 501 case(209) 502 if (.not.(spgroupma>=50 .and. 51>=spgroupma) ) spgrmatch=0 503 case(210) 504 if (.not.(spgroupma>=54 .and. 55>=spgroupma) ) spgrmatch=0 505 case(211) 506 if (.not.(spgroupma>=58 .and. 58>=spgroupma) ) spgrmatch=0 507 case(212) 508 if (.not.(spgroupma>=61 .and. 62>=spgroupma) ) spgrmatch=0 509 case(213) 510 if (.not.(spgroupma>=65 .and. 66>=spgroupma) ) spgrmatch=0 511 case(214) 512 if (.not.(spgroupma>=69 .and. 69>=spgroupma) ) spgrmatch=0 513 case(215) 514 if (.not.(spgroupma>=72 .and. 73>=spgroupma) ) spgrmatch=0 515 case(216) 516 if (.not.(spgroupma>=76 .and. 77>=spgroupma) ) spgrmatch=0 517 case(217) 518 if (.not.(spgroupma>=80 .and. 80>=spgroupma) ) spgrmatch=0 519 case(218) 520 if (.not.(spgroupma>=83 .and. 84>=spgroupma) ) spgrmatch=0 521 case(219) 522 if (.not.(spgroupma>=87 .and. 88>=spgroupma) ) spgrmatch=0 523 case(220) 524 if (.not.(spgroupma>=91 .and. 91>=spgroupma) ) spgrmatch=0 525 case(221) 526 if (.not.(spgroupma>=94 .and. 97>=spgroupma) ) spgrmatch=0 527 case(222) 528 if (.not.(spgroupma>=100 .and. 103>=spgroupma) ) spgrmatch=0 529 case(223) 530 if (.not.(spgroupma>=106 .and. 109>=spgroupma) ) spgrmatch=0 531 case(224) 532 if (.not.(spgroupma>=112 .and. 115>=spgroupma) ) spgrmatch=0 533 case(225) 534 if (.not.(spgroupma>=118 .and. 121>=spgroupma) ) spgrmatch=0 535 case(226) 536 if (.not.(spgroupma>=124 .and. 127>=spgroupma) ) spgrmatch=0 537 case(227) 538 if (.not.(spgroupma>=130 .and. 133>=spgroupma) ) spgrmatch=0 539 case(228) 540 if (.not.(spgroupma>=136 .and. 139>=spgroupma) ) spgrmatch=0 541 case(229) 542 if (.not.(spgroupma>=142 .and. 144>=spgroupma) ) spgrmatch=0 543 case(230) 544 if (.not.(spgroupma>=147 .and. 149>=spgroupma) ) spgrmatch=0 545 case default 546 write(message, '(3a,i8,4a)' )& 547 & 'The non-magnetic spacegroup is not specified ',ch10,& 548 & 'while the magnetic space group is specified, spgroupma= ',spgroupma,ch10,& 549 & 'This is not allowed. ',ch10,& 550 & 'Action: specify spgroup in the input file.' 551 MSG_ERROR(message) 552 end select 553 554 if (spgrmatch==0) then 555 write(message, '(a,i8,a,a,i8,4a)' )& 556 & 'mismatch between the non-magnetic spacegroup ',spgroup,ch10,& 557 & 'and the magnetic space group ',spgroupma,ch10,& 558 & 'This is not allowed. ',ch10,& 559 & 'Action: modify spgroup or spgroupma in the input file.' 560 MSG_ERROR(message) 561 end if 562 563 !DEBUG 564 !write(std_out,*) ' gensymshub, after check the spgroup ... ' 565 !ENDDEBUG 566 567 !Assign the magnetic Bravais lattice type from the magnetic space group number 568 !As the magnetic space group number begins with 1 for EACH crystal system 569 !we must first make our choice as a function of spgroup 570 571 !Convention : 572 !brvlttbw = input variable giving Bravais black-and-white translation 573 !(from 1 to 8 : Shubnikov type IV space group) 574 !1,2,3 = 1/2 translation along a, b, or c respectively 575 !4,5,6 = translation corresponding to A,B,C centering, respectively 576 !7 = I centering 577 !8 = (1/2 0 0) centering corresponding to normal F lattice 578 !9 -> Shubnikov type III space group 579 !Note that the use of the table 7.3 (p585) of Bradney and Cracknell 580 !for the definition of the translation vectors 581 !is extremely confusing, due to the strange choice of basis 582 !vectors of table 3.1. 583 !See table 7.4 (p588) for the spgroupma interpretation 584 585 select case(spgroup) 586 case(1,2) !Triclinic 587 select case(spgroupma) 588 case(6) 589 brvlttbw=9 !ShubIII 590 case(3,7) 591 brvlttbw=1 !Ps (note that it is not body centered, according to Table 592 ! 7.3 of Bradley and Cracknell) 593 end select 594 case(3:15) !Monoclinic 595 select case(spgroupma) 596 case(3,9,15,20,26,34,39,44:46,52:54,60:62,67:69,77:79,87:89) 597 brvlttbw=9 !ShubIII 598 case(4,10,17,21,27,36,41,47,55,64,70,80,91) 599 brvlttbw=1 !a 600 case(5,11,22,29,48,56,71,81) 601 brvlttbw=2 !b 602 case(16,28,35,40,63,72,82,90) 603 brvlttbw=3 !c 604 case(31,73,83) 605 brvlttbw=4 !A 606 case(6,12,23,30,49,57,74,84) 607 brvlttbw=6 !C 608 end select 609 case(16:74) !Orthorhombic 610 select case(spgroupma) 611 case(3,9,10,18,19,27,33,34,40,41,47,51,55,59,60,68:70,& 612 & 80,81,89:91,101:103,113:115,125:127,137,138,146:148,158,159,& 613 & 167,168,174:176,182,183,189:191,197:199,205:207,213:215,221,& 614 & 222,226,227,231,232,237,238,243:245,251:253,259:261,267:271,& 615 & 279:283,291:297,307:313,323:329,339:345,355:359,367:371,& 616 & 379:385,395:399,407:411,419:425,435:437,443:449,459:465,& 617 & 471:477,483:487,493:497,503:507,513:517,523:525,529:531,& 618 & 535:537,541:545,550:552,556:560) 619 brvlttbw=9 !ShubIII 620 case(4,11,20,28,36,43,62,71,83,92,104,116,128,140,& 621 & 149,160,170,178,185,192,200,208,216,234,240,247,254,262,272,& 622 & 284,298,314,330,346,360,372,386,400,412,426,438,450,467,& 623 & 479,489,499,509,519,547,562) 624 brvlttbw=1 !a 625 case(72,93,105,117,129,150,248,299,315,331,347,387,427,451) 626 brvlttbw=2 !b 627 case(12,21,35,42,52,56,61,73,82,94,106,118,130,139,& 628 & 151,161,169,177,184,193,201,209,217,233,239,246,273,285,300,& 629 & 316,332,348,361,373,388,401,413,428,452,466,478,488,498,& 630 & 508,518,538,546,553,561) 631 brvlttbw=3 !c 632 case(13,22,37,44,64,74,85,95,107,119,131,142,152,162,& 633 & 171,179,186,274,286,301,317,333,349,362,374,389,402,414,429,& 634 & 453,468,480,490,500,510,520) 635 brvlttbw=4 !A 636 case(75,96,108,120,132,153,302,318,334,350,390,430,454) 637 brvlttbw=5 !B 638 case(5,14,23,29,63,76,84,97,109,121,133,141,154,163,& 639 & 194,202,210,218,255,263,275,287,303,319,335,351,363,375,& 640 & 391,403,415,431,439,455) 641 brvlttbw=6 !C 642 case(6,15,24,30,65,77,86,98,110,122,134,143,155,164,256,& 643 & 264,276,288,304,320,336,352,364,376,392,404,416,432,440,456) 644 brvlttbw=7 !I 645 case(48,223,228,526,532) 646 brvlttbw=8 !s 647 end select 648 case(75:142) !Tetragonal 649 select case(spgroupma) 650 case(3,9,15,21,27,31,35,41,45:47,53:55,61:63,69:71,& 651 & 77:79,83:85,89:91,97:99,105:107,113:115,121:123,129:131,& 652 & 137:139,145:147,153:155,159:161,165:167,173:175,181:183,& 653 & 189:191,197:199,205:207,213:215,221:223,229:231,235:237,& 654 & 241:243,247:249,253:255,261:263,269:271,277:279,285:287,& 655 & 293:295,301:303,309:311,317:319,323:325,329:331,335:337,& 656 & 341:347,353:359,365:371,377:383,389:395,401:407,413:419,& 657 & 425:431,437:443,449:455,461:467,473:479,485:491,497:503,& 658 & 509:515,521:527,533:539,543:549,553:559,563:569) 659 brvlttbw=9 !ShubIII 660 case(4,10,16,22,28,32,36,42,48,56,64,72,80,86,92,100,& 661 & 108,116,124,132,140,148,156,162,168,176,184,192,200,208,& 662 & 216,224,232,238,244,250,256,264,272,280,288,296,304,312,& 663 & 320,326,332,338,348,360,372,384,396,408,420,432,444,456,& 664 & 468,480,492,504,516,528,540,550,560,570) 665 brvlttbw=3 !c 666 case(5,11,17,23,37,49,57,65,73,93,101,109,117,125,133,141,& 667 & 149,169,177,185,193,201,209,217,225,257,265,273,281,289,297,& 668 & 305,313,349,361,373,385,397,409,421,433,445,457,469,481,493,& 669 & 505,517,529) 670 brvlttbw=6 !C 671 case(6,12,18,24,38,50,58,66,74,94,102,110,118,126,134,142,& 672 & 150,170,178,186,194,202,210,218,226,258,266,274,282,290,298,& 673 & 306,314,350,362,374,386,398,410,422,434,446,458,470,482,494,& 674 & 506,518,530) 675 brvlttbw=7 !I 676 end select 677 case(143:194) !Hexagonal 678 select case(spgroupma) 679 case(15,19,23,27,31,35,39,43,47,51,55,59,63,67,71,& 680 & 75:77,81:83,87:89,93:95,99:101,105:107,111,115,119,123,127,& 681 & 131,135,139:141,145:147,151:153,157:159,163:165,169:171,& 682 & 175:177,181:183,187:189,193:195,199:201,205:207,211:213,& 683 & 217:219,223:225,229:231,235:241,245:251,255:261,265:271) 684 brvlttbw=9 !ShubIII 685 case(3,6,9,16,24,28,32,36,40,44,52,56,60,64,78,84,& 686 & 90,96,112,116,120,124,128,132,136,142,148,154,160,166,172,& 687 & 178,184,190,196,202,208,214,220,226,232,242,252,262,272) 688 brvlttbw=6 !C 689 case(12,20,48,68,72,102,108) 690 brvlttbw=7 !I 691 end select 692 case(195:230) !Cubic 693 select case(spgroupma) 694 case(16,20,24,28,32,35,39,42,46,50,54,58,61,65,69,& 695 & 72,76,80,83,87,91,94:96,100:102,106:108,112:114,118:120,& 696 & 124:126,130:132,136:138,142:144,147:149) 697 brvlttbw=9 !ShubIII 698 case(3,11,17,21,36,43,47,62,66,73,84,97,103,109,115) 699 brvlttbw=7 !I 700 case(6,25,29,51,55,77,88,121,127,133,139) 701 brvlttbw=8 !s 702 end select 703 end select 704 705 if(brvlttbw==9)shubnikov=3 !Shubnikov type III 706 if(brvlttbw>=1 .and. brvlttbw<=8) shubnikov=4 !Shubnikov type IV 707 708 genafm(:)=zero 709 if(shubnikov==4)then 710 if(brvlttbw==1)genafm(1)=half 711 if(brvlttbw==2)genafm(2)=half 712 if(brvlttbw==3)genafm(3)=half 713 if(brvlttbw>=4 .and. brvlttbw<=8)then 714 genafm(:)=half 715 if(brvlttbw==4)genafm(1)=zero 716 if(brvlttbw==5)genafm(2)=zero 717 if(brvlttbw==6)genafm(3)=zero 718 end if 719 end if 720 721 !DEBUG 722 !write(std_out,*) 'gensymshub : end ' 723 !write(std_out,*) 'gensymshub, shubnikov =',shubnikov 724 !ENDDEBUG 725 end subroutine gensymshub