TABLE OF CONTENTS


ABINIT/m_spgbuilder [ Modules ]

[ Top ] [ Modules ]

NAME

 m_spgbuilder

FUNCTION

  Spacegroup builder.

COPYRIGHT

  Copyright (C) 2008-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 .

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_spgbuilder
28 
29  use defs_basis
30  use m_errors
31  use m_abicore
32 
33  use m_symtk,   only : chkgrp, print_symmetries
34  use m_symsg,   only : symsgcube, symsghexa, symsgmono, symsgortho, symsgtetra
35 
36  implicit none
37 
38  private

m_spgbuilder/gensymshub [ Functions ]

[ Top ] [ m_spgbuilder ] [ 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).

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

 607 subroutine gensymshub(genafm,spgroup,spgroupma,shubnikov)
 608 
 609 
 610 !This section has been created automatically by the script Abilint (TD).
 611 !Do not modify the following lines by hand.
 612 #undef ABI_FUNC
 613 #define ABI_FUNC 'gensymshub'
 614 !End of the abilint section
 615 
 616  implicit none
 617 
 618 !Arguments ------------------------------------
 619 !scalars
 620  integer,intent(in) :: spgroup,spgroupma
 621  integer,intent(out) :: shubnikov
 622 !arrays
 623  real(dp),intent(out) :: genafm(3)
 624 
 625 !Local variables ------------------------------
 626 !scalars
 627  integer :: brvlttbw=0,spgrmatch=1
 628  character(len=500) :: message
 629 
 630 ! *************************************************************************
 631 
 632 !List of the input parameters
 633 !DEBUG
 634 !write(std_out,*)
 635 !write(std_out,*)' gensymshub : enter with:'
 636 !write(std_out,*)' brvlttbw  = ',brvlttbw
 637 !write(std_out,*)' spgroup = ',spgroup
 638 !write(std_out,*)' spgroupma = ',spgroupma
 639 !ENDDEBUG
 640 
 641 !Test for consistency the magnetic and non-magnetic space group
 642  select case (spgroup)
 643  case(1)
 644    if (.not.(spgroupma>=3  .and.  3 >=spgroupma) ) spgrmatch=0
 645  case(2)
 646    if (.not.(spgroupma>=6  .and.  7 >=spgroupma) ) spgrmatch=0
 647  case(3)
 648    if (.not.(spgroupma>=3   .and.  6 >=spgroupma) ) spgrmatch=0
 649  case(4)
 650    if (.not.(spgroupma>=9   .and.  12 >=spgroupma) ) spgrmatch=0
 651  case(5)
 652    if (.not.(spgroupma>=15   .and.  17 >=spgroupma) ) spgrmatch=0
 653  case(6)
 654    if (.not.(spgroupma>=20   .and.  23 >=spgroupma) ) spgrmatch=0
 655  case(7)
 656    if (.not.(spgroupma>=26   .and.  31 >=spgroupma) ) spgrmatch=0
 657  case(8)
 658    if (.not.(spgroupma>=34   .and.  36 >=spgroupma) ) spgrmatch=0
 659  case(9)
 660    if (.not.(spgroupma>=39   .and.  41 >=spgroupma) ) spgrmatch=0
 661  case(10)
 662    if (.not.(spgroupma>=44   .and.  49 >=spgroupma) ) spgrmatch=0
 663  case(11)
 664    if (.not.(spgroupma>=52   .and.  57 >=spgroupma) ) spgrmatch=0
 665  case(12)
 666    if (.not.(spgroupma>=60   .and.  64 >=spgroupma) ) spgrmatch=0
 667  case(13)
 668    if (.not.(spgroupma>=67   .and.  74 >=spgroupma) ) spgrmatch=0
 669  case(14)
 670    if (.not.(spgroupma>=77   .and.  84 >=spgroupma) ) spgrmatch=0
 671  case(15)
 672    if (.not.(spgroupma>=87   .and.  91 >=spgroupma) ) spgrmatch=0
 673  case(16)
 674    if (.not.(spgroupma>= 3  .and.   6>=spgroupma) ) spgrmatch=0
 675  case(17)
 676    if (.not.(spgroupma>= 9  .and.  15 >=spgroupma) ) spgrmatch=0
 677  case(18)
 678    if (.not.(spgroupma>=18   .and.  24 >=spgroupma) ) spgrmatch=0
 679  case(19)
 680    if (.not.(spgroupma>=27   .and.  30 >=spgroupma) ) spgrmatch=0
 681  case(20)
 682    if (.not.(spgroupma>=33   .and.  37 >=spgroupma) ) spgrmatch=0
 683  case(21)
 684    if (.not.(spgroupma>=40   .and.  44 >=spgroupma) ) spgrmatch=0
 685  case(22)
 686    if (.not.(spgroupma>=47   .and.  48 >=spgroupma) ) spgrmatch=0
 687  case(23)
 688    if (.not.(spgroupma>=51   .and.  52 >=spgroupma) ) spgrmatch=0
 689  case(24)
 690    if (.not.(spgroupma>=55   .and.  56 >=spgroupma) ) spgrmatch=0
 691  case(25)
 692    if (.not.(spgroupma>=59   .and.  65 >=spgroupma) ) spgrmatch=0
 693  case(26)
 694    if (.not.(spgroupma>=68   .and.  77 >=spgroupma) ) spgrmatch=0
 695  case(27)
 696    if (.not.(spgroupma>=80   .and.  86 >=spgroupma) ) spgrmatch=0
 697  case(28)
 698    if (.not.(spgroupma>=89   .and.  98 >=spgroupma) ) spgrmatch=0
 699  case(29)
 700    if (.not.(spgroupma>=101  .and.  110 >=spgroupma) ) spgrmatch=0
 701  case(30)
 702    if (.not.(spgroupma>=113   .and. 122  >=spgroupma) ) spgrmatch=0
 703  case(31)
 704    if (.not.(spgroupma>=125   .and. 134  >=spgroupma) ) spgrmatch=0
 705  case(32)
 706    if (.not.(spgroupma>=137   .and. 143  >=spgroupma) ) spgrmatch=0
 707  case(33)
 708    if (.not.(spgroupma>=146   .and.  155 >=spgroupma) ) spgrmatch=0
 709  case(34)
 710    if (.not.(spgroupma>=158   .and.  164 >=spgroupma) ) spgrmatch=0
 711  case(35)
 712    if (.not.(spgroupma>=167   .and.   171>=spgroupma) ) spgrmatch=0
 713  case(36)
 714    if (.not.(spgroupma>=174   .and.   179>=spgroupma) ) spgrmatch=0
 715  case(37)
 716    if (.not.(spgroupma>=182   .and.  186 >=spgroupma) ) spgrmatch=0
 717  case(38)
 718    if (.not.(spgroupma>=189   .and.  194 >=spgroupma) ) spgrmatch=0
 719  case(39)
 720    if (.not.(spgroupma>=197   .and.  202 >=spgroupma) ) spgrmatch=0
 721  case(40)
 722    if (.not.(spgroupma>=205   .and.  210 >=spgroupma) ) spgrmatch=0
 723  case(41)
 724    if (.not.(spgroupma>=213   .and.  218 >=spgroupma) ) spgrmatch=0
 725  case(42)
 726    if (.not.(spgroupma>=221   .and.   223>=spgroupma) ) spgrmatch=0
 727  case(43)
 728    if (.not.(spgroupma>=226   .and.   228>=spgroupma) ) spgrmatch=0
 729  case(44)
 730    if (.not.(spgroupma>=231   .and.  234 >=spgroupma) ) spgrmatch=0
 731  case(45)
 732    if (.not.(spgroupma>=237   .and.  240 >=spgroupma) ) spgrmatch=0
 733  case(46)
 734    if (.not.(spgroupma>=243   .and.  248 >=spgroupma) ) spgrmatch=0
 735  case(47)
 736    if (.not.(spgroupma>=251   .and.  256 >=spgroupma) ) spgrmatch=0
 737  case(48)
 738    if (.not.(spgroupma>=259   .and.  264 >=spgroupma) ) spgrmatch=0
 739  case(49)
 740    if (.not.(spgroupma>=267   .and.  276 >=spgroupma) ) spgrmatch=0
 741  case(50)
 742    if (.not.(spgroupma>=279   .and.  288 >=spgroupma) ) spgrmatch=0
 743  case(51)
 744    if (.not.(spgroupma>=291   .and.  304 >=spgroupma) ) spgrmatch=0
 745  case(52)
 746    if (.not.(spgroupma>=307   .and.  320 >=spgroupma) ) spgrmatch=0
 747  case(53)
 748    if (.not.(spgroupma>=323   .and.  336 >=spgroupma) ) spgrmatch=0
 749  case(54)
 750    if (.not.(spgroupma>=339   .and.  352 >=spgroupma) ) spgrmatch=0
 751  case(55)
 752    if (.not.(spgroupma>=355   .and.  364 >=spgroupma) ) spgrmatch=0
 753  case(56)
 754    if (.not.(spgroupma>=367   .and.  376 >=spgroupma) ) spgrmatch=0
 755  case(57)
 756    if (.not.(spgroupma>=379   .and.  392 >=spgroupma) ) spgrmatch=0
 757  case(58)
 758    if (.not.(spgroupma>=395   .and.  404 >=spgroupma) ) spgrmatch=0
 759  case(59)
 760    if (.not.(spgroupma>=407   .and.  416 >=spgroupma) ) spgrmatch=0
 761  case(60)
 762    if (.not.(spgroupma>=419   .and.  432 >=spgroupma) ) spgrmatch=0
 763  case(61)
 764    if (.not.(spgroupma>=435   .and.  440 >=spgroupma) ) spgrmatch=0
 765  case(62)
 766    if (.not.(spgroupma>=443   .and.  456 >=spgroupma) ) spgrmatch=0
 767  case(63)
 768    if (.not.(spgroupma>=459   .and.  468 >=spgroupma) ) spgrmatch=0
 769  case(64)
 770    if (.not.(spgroupma>=471   .and.  480 >=spgroupma) ) spgrmatch=0
 771  case(65)
 772    if (.not.(spgroupma>=483   .and.  490 >=spgroupma) ) spgrmatch=0
 773  case(66)
 774    if (.not.(spgroupma>=493   .and.  500 >=spgroupma) ) spgrmatch=0
 775  case(67)
 776    if (.not.(spgroupma>=503   .and.  510 >=spgroupma) ) spgrmatch=0
 777  case(68)
 778    if (.not.(spgroupma>=513   .and.  520 >=spgroupma) ) spgrmatch=0
 779  case(69)
 780    if (.not.(spgroupma>=523   .and.  526 >=spgroupma) ) spgrmatch=0
 781  case(70)
 782    if (.not.(spgroupma>=529   .and.  532 >=spgroupma) ) spgrmatch=0
 783  case(71)
 784    if (.not.(spgroupma>=535   .and.  538 >=spgroupma) ) spgrmatch=0
 785  case(72)
 786    if (.not.(spgroupma>=541   .and.  547 >=spgroupma) ) spgrmatch=0
 787  case(73)
 788    if (.not.(spgroupma>=550   .and.  553 >=spgroupma) ) spgrmatch=0
 789  case(74)
 790    if (.not.(spgroupma>=556   .and.  562 >=spgroupma) ) spgrmatch=0
 791  case(75)
 792    if (.not.(spgroupma>= 3  .and.   6>=spgroupma) ) spgrmatch=0
 793  case(76)
 794    if (.not.(spgroupma>= 9  .and.  12 >=spgroupma) ) spgrmatch=0
 795  case(77)
 796    if (.not.(spgroupma>= 15  .and.   18>=spgroupma) ) spgrmatch=0
 797  case(78)
 798    if (.not.(spgroupma>= 21  .and.   24>=spgroupma) ) spgrmatch=0
 799  case(79)
 800    if (.not.(spgroupma>= 27  .and.  28 >=spgroupma) ) spgrmatch=0
 801  case(80)
 802    if (.not.(spgroupma>= 31  .and.  32 >=spgroupma) ) spgrmatch=0
 803  case(81)
 804    if (.not.(spgroupma>= 35  .and.  38 >=spgroupma) ) spgrmatch=0
 805  case(82)
 806    if (.not.(spgroupma>= 41  .and.  42 >=spgroupma) ) spgrmatch=0
 807  case(83)
 808    if (.not.(spgroupma>=45   .and.   50>=spgroupma) ) spgrmatch=0
 809  case(84)
 810    if (.not.(spgroupma>=53   .and.   58>=spgroupma) ) spgrmatch=0
 811  case(85)
 812    if (.not.(spgroupma>=61   .and.   66>=spgroupma) ) spgrmatch=0
 813  case(86)
 814    if (.not.(spgroupma>=69   .and.   74>=spgroupma) ) spgrmatch=0
 815  case(87)
 816    if (.not.(spgroupma>=77   .and.   80>=spgroupma) ) spgrmatch=0
 817  case(88)
 818    if (.not.(spgroupma>=83   .and.   86>=spgroupma) ) spgrmatch=0
 819  case(89)
 820    if (.not.(spgroupma>=89   .and.   94>=spgroupma) ) spgrmatch=0
 821  case(90)
 822    if (.not.(spgroupma>=97   .and.   102>=spgroupma) ) spgrmatch=0
 823  case(91)
 824    if (.not.(spgroupma>=105   .and.   110>=spgroupma) ) spgrmatch=0
 825  case(92)
 826    if (.not.(spgroupma>=113   .and.   118>=spgroupma) ) spgrmatch=0
 827  case(93)
 828    if (.not.(spgroupma>=121   .and.   126>=spgroupma) ) spgrmatch=0
 829  case(94)
 830    if (.not.(spgroupma>=129   .and.   134>=spgroupma) ) spgrmatch=0
 831  case(95)
 832    if (.not.(spgroupma>=137   .and.   142>=spgroupma) ) spgrmatch=0
 833  case(96)
 834    if (.not.(spgroupma>=145   .and.   150>=spgroupma) ) spgrmatch=0
 835  case(97)
 836    if (.not.(spgroupma>=153   .and.   156>=spgroupma) ) spgrmatch=0
 837  case(98)
 838    if (.not.(spgroupma>=159   .and.   162>=spgroupma) ) spgrmatch=0
 839  case(99)
 840    if (.not.(spgroupma>=165   .and.   170>=spgroupma) ) spgrmatch=0
 841  case(100)
 842    if (.not.(spgroupma>=173   .and.   178>=spgroupma) ) spgrmatch=0
 843  case(101)
 844    if (.not.(spgroupma>=181   .and.   186>=spgroupma) ) spgrmatch=0
 845  case(102)
 846    if (.not.(spgroupma>=189   .and.   194>=spgroupma) ) spgrmatch=0
 847  case(103)
 848    if (.not.(spgroupma>=197   .and.   202>=spgroupma) ) spgrmatch=0
 849  case(104)
 850    if (.not.(spgroupma>=205   .and.   210>=spgroupma) ) spgrmatch=0
 851  case(105)
 852    if (.not.(spgroupma>=213   .and.   218>=spgroupma) ) spgrmatch=0
 853  case(106)
 854    if (.not.(spgroupma>=221   .and.   226>=spgroupma) ) spgrmatch=0
 855  case(107)
 856    if (.not.(spgroupma>=229  .and.   232>=spgroupma) ) spgrmatch=0
 857  case(108)
 858    if (.not.(spgroupma>=235   .and.   238>=spgroupma) ) spgrmatch=0
 859  case(109)
 860    if (.not.(spgroupma>=241   .and.   244>=spgroupma) ) spgrmatch=0
 861  case(110)
 862    if (.not.(spgroupma>=247   .and.   250>=spgroupma) ) spgrmatch=0
 863  case(111)
 864    if (.not.(spgroupma>=253   .and.   258>=spgroupma) ) spgrmatch=0
 865  case(112)
 866    if (.not.(spgroupma>=261   .and.   266>=spgroupma) ) spgrmatch=0
 867  case(113)
 868    if (.not.(spgroupma>=269   .and.   274>=spgroupma) ) spgrmatch=0
 869  case(114)
 870    if (.not.(spgroupma>=277   .and.   282>=spgroupma) ) spgrmatch=0
 871  case(115)
 872    if (.not.(spgroupma>=285   .and.   290>=spgroupma) ) spgrmatch=0
 873  case(116)
 874    if (.not.(spgroupma>=293   .and.   298>=spgroupma) ) spgrmatch=0
 875  case(117)
 876    if (.not.(spgroupma>=301   .and.   306>=spgroupma) ) spgrmatch=0
 877  case(118)
 878    if (.not.(spgroupma>=309   .and.   314>=spgroupma) ) spgrmatch=0
 879  case(119)
 880    if (.not.(spgroupma>=317   .and.   320>=spgroupma) ) spgrmatch=0
 881  case(120)
 882    if (.not.(spgroupma>=323   .and.   326>=spgroupma) ) spgrmatch=0
 883  case(121)
 884    if (.not.(spgroupma>=329   .and.   332>=spgroupma) ) spgrmatch=0
 885  case(122)
 886    if (.not.(spgroupma>=335   .and.   338>=spgroupma) ) spgrmatch=0
 887  case(123)
 888    if (.not.(spgroupma>=341   .and.   350>=spgroupma) ) spgrmatch=0
 889  case(124)
 890    if (.not.(spgroupma>=353   .and.   362>=spgroupma) ) spgrmatch=0
 891  case(125)
 892    if (.not.(spgroupma>=365   .and.   374>=spgroupma) ) spgrmatch=0
 893  case(126)
 894    if (.not.(spgroupma>=377   .and.   386>=spgroupma) ) spgrmatch=0
 895  case(127)
 896    if (.not.(spgroupma>=389   .and.   398>=spgroupma) ) spgrmatch=0
 897  case(128)
 898    if (.not.(spgroupma>=401   .and.   410>=spgroupma) ) spgrmatch=0
 899  case(129)
 900    if (.not.(spgroupma>=413   .and.   422>=spgroupma) ) spgrmatch=0
 901  case(130)
 902    if (.not.(spgroupma>=425   .and.   434>=spgroupma) ) spgrmatch=0
 903  case(131)
 904    if (.not.(spgroupma>=437   .and.   446>=spgroupma) ) spgrmatch=0
 905  case(132)
 906    if (.not.(spgroupma>=449   .and.   458>=spgroupma) ) spgrmatch=0
 907  case(133)
 908    if (.not.(spgroupma>=461   .and.   470>=spgroupma) ) spgrmatch=0
 909  case(134)
 910    if (.not.(spgroupma>=473   .and.   482>=spgroupma) ) spgrmatch=0
 911  case(135)
 912    if (.not.(spgroupma>=485   .and.   494>=spgroupma) ) spgrmatch=0
 913  case(136)
 914    if (.not.(spgroupma>=497  .and.   506>=spgroupma) ) spgrmatch=0
 915  case(137)
 916    if (.not.(spgroupma>=509   .and.   518>=spgroupma) ) spgrmatch=0
 917  case(138)
 918    if (.not.(spgroupma>=521   .and.   530>=spgroupma) ) spgrmatch=0
 919  case(139)
 920    if (.not.(spgroupma>=533   .and.   540>=spgroupma) ) spgrmatch=0
 921  case(140)
 922    if (.not.(spgroupma>=543   .and.   550>=spgroupma) ) spgrmatch=0
 923  case(141)
 924    if (.not.(spgroupma>=553   .and.   560>=spgroupma) ) spgrmatch=0
 925  case(142)
 926    if (.not.(spgroupma>=563   .and.   570>=spgroupma) ) spgrmatch=0
 927  case(143)
 928    if (.not.(spgroupma>=3   .and.   3>=spgroupma) ) spgrmatch=0
 929  case(144)
 930    if (.not.(spgroupma>=6   .and.   6>=spgroupma) ) spgrmatch=0
 931  case(145)
 932    if (.not.(spgroupma>=9   .and.   9>=spgroupma) ) spgrmatch=0
 933  case(146)
 934    if (.not.(spgroupma>=12   .and.   12>=spgroupma) ) spgrmatch=0
 935  case(147)
 936    if (.not.(spgroupma>=15   .and.   16>=spgroupma) ) spgrmatch=0
 937  case(148)
 938    if (.not.(spgroupma>=19   .and.   20>=spgroupma) ) spgrmatch=0
 939  case(149)
 940    if (.not.(spgroupma>=23   .and.   24>=spgroupma) ) spgrmatch=0
 941  case(150)
 942    if (.not.(spgroupma>=27   .and.   28>=spgroupma) ) spgrmatch=0
 943  case(151)
 944    if (.not.(spgroupma>=31   .and.   32>=spgroupma) ) spgrmatch=0
 945  case(152)
 946    if (.not.(spgroupma>=35   .and.   36>=spgroupma) ) spgrmatch=0
 947  case(153)
 948    if (.not.(spgroupma>=39   .and.   40>=spgroupma) ) spgrmatch=0
 949  case(154)
 950    if (.not.(spgroupma>=43   .and.   44>=spgroupma) ) spgrmatch=0
 951  case(155)
 952    if (.not.(spgroupma>=47   .and.   48>=spgroupma) ) spgrmatch=0
 953  case(156)
 954    if (.not.(spgroupma>=51   .and.   52>=spgroupma) ) spgrmatch=0
 955  case(157)
 956    if (.not.(spgroupma>=55   .and.   56>=spgroupma) ) spgrmatch=0
 957  case(158)
 958    if (.not.(spgroupma>=59   .and.   60>=spgroupma) ) spgrmatch=0
 959  case(159)
 960    if (.not.(spgroupma>=63   .and.   64>=spgroupma) ) spgrmatch=0
 961  case(160)
 962    if (.not.(spgroupma>=67   .and.   68>=spgroupma) ) spgrmatch=0
 963  case(161)
 964    if (.not.(spgroupma>=71   .and.   72>=spgroupma) ) spgrmatch=0
 965  case(162)
 966    if (.not.(spgroupma>=75   .and.   78>=spgroupma) ) spgrmatch=0
 967  case(163)
 968    if (.not.(spgroupma>=81   .and.   84>=spgroupma) ) spgrmatch=0
 969  case(164)
 970    if (.not.(spgroupma>=87   .and.   90>=spgroupma) ) spgrmatch=0
 971  case(165)
 972    if (.not.(spgroupma>=93   .and.   96>=spgroupma) ) spgrmatch=0
 973  case(166)
 974    if (.not.(spgroupma>=99   .and.   102>=spgroupma) ) spgrmatch=0
 975  case(167)
 976    if (.not.(spgroupma>=105   .and.   108>=spgroupma) ) spgrmatch=0
 977  case(168)
 978    if (.not.(spgroupma>=111   .and.   112>=spgroupma) ) spgrmatch=0
 979  case(169)
 980    if (.not.(spgroupma>=115   .and.   116>=spgroupma) ) spgrmatch=0
 981  case(170)
 982    if (.not.(spgroupma>=119   .and.   120>=spgroupma) ) spgrmatch=0
 983  case(171)
 984    if (.not.(spgroupma>=123   .and.   124>=spgroupma) ) spgrmatch=0
 985  case(172)
 986    if (.not.(spgroupma>=127   .and.   128>=spgroupma) ) spgrmatch=0
 987  case(173)
 988    if (.not.(spgroupma>=131   .and.   132>=spgroupma) ) spgrmatch=0
 989  case(174)
 990    if (.not.(spgroupma>=135   .and.   136>=spgroupma) ) spgrmatch=0
 991  case(175)
 992    if (.not.(spgroupma>=139   .and.   142>=spgroupma) ) spgrmatch=0
 993  case(176)
 994    if (.not.(spgroupma>=145   .and.   148>=spgroupma) ) spgrmatch=0
 995  case(177)
 996    if (.not.(spgroupma>=151   .and.   154>=spgroupma) ) spgrmatch=0
 997  case(178)
 998    if (.not.(spgroupma>=157   .and.   160>=spgroupma) ) spgrmatch=0
 999  case(179)
1000    if (.not.(spgroupma>=163   .and.   166>=spgroupma) ) spgrmatch=0
1001  case(180)
1002    if (.not.(spgroupma>=169   .and.   172>=spgroupma) ) spgrmatch=0
1003  case(181)
1004    if (.not.(spgroupma>=175   .and.   178>=spgroupma) ) spgrmatch=0
1005  case(182)
1006    if (.not.(spgroupma>=181   .and.   184>=spgroupma) ) spgrmatch=0
1007  case(183)
1008    if (.not.(spgroupma>=187   .and.   190>=spgroupma) ) spgrmatch=0
1009  case(184)
1010    if (.not.(spgroupma>=193   .and.   196>=spgroupma) ) spgrmatch=0
1011  case(185)
1012    if (.not.(spgroupma>=199   .and.   202>=spgroupma) ) spgrmatch=0
1013  case(186)
1014    if (.not.(spgroupma>=205   .and.   208>=spgroupma) ) spgrmatch=0
1015  case(187)
1016    if (.not.(spgroupma>=211   .and.   214>=spgroupma) ) spgrmatch=0
1017  case(188)
1018    if (.not.(spgroupma>=217   .and.   220>=spgroupma) ) spgrmatch=0
1019  case(189)
1020    if (.not.(spgroupma>=223   .and.   226>=spgroupma) ) spgrmatch=0
1021  case(190)
1022    if (.not.(spgroupma>=229   .and.   232>=spgroupma) ) spgrmatch=0
1023  case(191)
1024    if (.not.(spgroupma>=235   .and.   242>=spgroupma) ) spgrmatch=0
1025  case(192)
1026    if (.not.(spgroupma>=245   .and.   252>=spgroupma) ) spgrmatch=0
1027  case(193)
1028    if (.not.(spgroupma>=255   .and.   262>=spgroupma) ) spgrmatch=0
1029  case(194)
1030    if (.not.(spgroupma>=265   .and.   272>=spgroupma) ) spgrmatch=0
1031  case(195)
1032    if (.not.(spgroupma>=3   .and.   3>=spgroupma) ) spgrmatch=0
1033  case(196)
1034    if (.not.(spgroupma>=6   .and.   6>=spgroupma) ) spgrmatch=0
1035  case(197)
1036    spgrmatch=0
1037  case(198)
1038    if (.not.(spgroupma>=11  .and.   11>=spgroupma) ) spgrmatch=0
1039  case(199)
1040    spgrmatch=0
1041  case(200)
1042    if (.not.(spgroupma>=16   .and.   17>=spgroupma) ) spgrmatch=0
1043  case(201)
1044    if (.not.(spgroupma>=20   .and.   21>=spgroupma) ) spgrmatch=0
1045  case(202)
1046    if (.not.(spgroupma>=24   .and.   25>=spgroupma) ) spgrmatch=0
1047  case(203)
1048    if (.not.(spgroupma>=28   .and.   29>=spgroupma) ) spgrmatch=0
1049  case(204)
1050    if (.not.(spgroupma>=32   .and.   32>=spgroupma) ) spgrmatch=0
1051  case(205)
1052    if (.not.(spgroupma>=35   .and.   36>=spgroupma) ) spgrmatch=0
1053  case(206)
1054    if (.not.(spgroupma>=39   .and.   39>=spgroupma) ) spgrmatch=0
1055  case(207)
1056    if (.not.(spgroupma>=42   .and.   43>=spgroupma) ) spgrmatch=0
1057  case(208)
1058    if (.not.(spgroupma>=46   .and.   47>=spgroupma) ) spgrmatch=0
1059  case(209)
1060    if (.not.(spgroupma>=50   .and.   51>=spgroupma) ) spgrmatch=0
1061  case(210)
1062    if (.not.(spgroupma>=54   .and.   55>=spgroupma) ) spgrmatch=0
1063  case(211)
1064    if (.not.(spgroupma>=58   .and.   58>=spgroupma) ) spgrmatch=0
1065  case(212)
1066    if (.not.(spgroupma>=61  .and.   62>=spgroupma) ) spgrmatch=0
1067  case(213)
1068    if (.not.(spgroupma>=65   .and.   66>=spgroupma) ) spgrmatch=0
1069  case(214)
1070    if (.not.(spgroupma>=69   .and.   69>=spgroupma) ) spgrmatch=0
1071  case(215)
1072    if (.not.(spgroupma>=72   .and.   73>=spgroupma) ) spgrmatch=0
1073  case(216)
1074    if (.not.(spgroupma>=76   .and.   77>=spgroupma) ) spgrmatch=0
1075  case(217)
1076    if (.not.(spgroupma>=80   .and.   80>=spgroupma) ) spgrmatch=0
1077  case(218)
1078    if (.not.(spgroupma>=83   .and.   84>=spgroupma) ) spgrmatch=0
1079  case(219)
1080    if (.not.(spgroupma>=87   .and.   88>=spgroupma) ) spgrmatch=0
1081  case(220)
1082    if (.not.(spgroupma>=91   .and.   91>=spgroupma) ) spgrmatch=0
1083  case(221)
1084    if (.not.(spgroupma>=94   .and.   97>=spgroupma) ) spgrmatch=0
1085  case(222)
1086    if (.not.(spgroupma>=100  .and.   103>=spgroupma) ) spgrmatch=0
1087  case(223)
1088    if (.not.(spgroupma>=106   .and.   109>=spgroupma) ) spgrmatch=0
1089  case(224)
1090    if (.not.(spgroupma>=112   .and.   115>=spgroupma) ) spgrmatch=0
1091  case(225)
1092    if (.not.(spgroupma>=118   .and.   121>=spgroupma) ) spgrmatch=0
1093  case(226)
1094    if (.not.(spgroupma>=124   .and.   127>=spgroupma) ) spgrmatch=0
1095  case(227)
1096    if (.not.(spgroupma>=130   .and.   133>=spgroupma) ) spgrmatch=0
1097  case(228)
1098    if (.not.(spgroupma>=136   .and.   139>=spgroupma) ) spgrmatch=0
1099  case(229)
1100    if (.not.(spgroupma>=142   .and.   144>=spgroupma) ) spgrmatch=0
1101  case(230)
1102    if (.not.(spgroupma>=147   .and.   149>=spgroupma) ) spgrmatch=0
1103  case default
1104    write(message, '(3a,i8,4a)' )&
1105 &   'The non-magnetic spacegroup is not specified ',ch10,&
1106 &   'while the magnetic space group is specified, spgroupma= ',spgroupma,ch10,&
1107 &   'This is not allowed.  ',ch10,&
1108 &   'Action: specify spgroup in the input file.'
1109    MSG_ERROR(message)
1110  end select
1111 
1112  if (spgrmatch==0) then
1113    write(message, '(a,i8,a,a,i8,4a)' )&
1114 &   'mismatch between the non-magnetic spacegroup ',spgroup,ch10,&
1115 &   'and the magnetic space group ',spgroupma,ch10,&
1116 &   'This is not allowed.  ',ch10,&
1117 &   'Action: modify spgroup or spgroupma in the input file.'
1118    MSG_ERROR(message)
1119  end if
1120 
1121 !DEBUG
1122 !write(std_out,*) ' gensymshub, after check the spgroup ... '
1123 !ENDDEBUG
1124 
1125 !Assign the magnetic Bravais lattice type from the magnetic space group number
1126 !As the magnetic space group number begins with 1 for EACH crystal system
1127 !we must first make our choice as a function of spgroup
1128 
1129 !Convention :
1130 !brvlttbw = input variable giving Bravais black-and-white translation
1131 !(from 1 to 8 : Shubnikov type IV space group)
1132 !1,2,3 = 1/2 translation along a, b, or c respectively
1133 !4,5,6 = translation corresponding to A,B,C centering, respectively
1134 !7 = I centering
1135 !8 = (1/2 0 0) centering corresponding to normal F lattice
1136 !9 -> Shubnikov type III space group
1137 !Note that the use of the table 7.3 (p585) of Bradney and Cracknell
1138 !for the definition of the translation vectors
1139 !is extremely confusing, due to the strange choice of basis
1140 !vectors of table 3.1.
1141 !See table 7.4 (p588) for the spgroupma interpretation
1142 
1143  select case(spgroup)
1144  case(1,2)     !Triclinic
1145    select case(spgroupma)
1146    case(6)
1147      brvlttbw=9      !ShubIII
1148    case(3,7)
1149      brvlttbw=1      !Ps (note that it is not body centered, according to Table
1150 !        7.3 of Bradley and Cracknell)
1151    end select
1152  case(3:15)     !Monoclinic
1153    select case(spgroupma)
1154    case(3,9,15,20,26,34,39,44:46,52:54,60:62,67:69,77:79,87:89)
1155      brvlttbw=9      !ShubIII
1156    case(4,10,17,21,27,36,41,47,55,64,70,80,91)
1157      brvlttbw=1     !a
1158    case(5,11,22,29,48,56,71,81)
1159      brvlttbw=2     !b
1160    case(16,28,35,40,63,72,82,90)
1161      brvlttbw=3     !c
1162    case(31,73,83)
1163      brvlttbw=4     !A
1164    case(6,12,23,30,49,57,74,84)
1165      brvlttbw=6     !C
1166    end select
1167  case(16:74)     !Orthorhombic
1168    select case(spgroupma)
1169    case(3,9,10,18,19,27,33,34,40,41,47,51,55,59,60,68:70,&
1170 &     80,81,89:91,101:103,113:115,125:127,137,138,146:148,158,159,&
1171 &     167,168,174:176,182,183,189:191,197:199,205:207,213:215,221,&
1172 &     222,226,227,231,232,237,238,243:245,251:253,259:261,267:271,&
1173 &     279:283,291:297,307:313,323:329,339:345,355:359,367:371,&
1174 &     379:385,395:399,407:411,419:425,435:437,443:449,459:465,&
1175 &     471:477,483:487,493:497,503:507,513:517,523:525,529:531,&
1176 &     535:537,541:545,550:552,556:560)
1177      brvlttbw=9      !ShubIII
1178    case(4,11,20,28,36,43,62,71,83,92,104,116,128,140,&
1179 &     149,160,170,178,185,192,200,208,216,234,240,247,254,262,272,&
1180 &     284,298,314,330,346,360,372,386,400,412,426,438,450,467,&
1181 &     479,489,499,509,519,547,562)
1182      brvlttbw=1     !a
1183    case(72,93,105,117,129,150,248,299,315,331,347,387,427,451)
1184      brvlttbw=2     !b
1185    case(12,21,35,42,52,56,61,73,82,94,106,118,130,139,&
1186 &     151,161,169,177,184,193,201,209,217,233,239,246,273,285,300,&
1187 &     316,332,348,361,373,388,401,413,428,452,466,478,488,498,&
1188 &     508,518,538,546,553,561)
1189      brvlttbw=3     !c
1190    case(13,22,37,44,64,74,85,95,107,119,131,142,152,162,&
1191 &     171,179,186,274,286,301,317,333,349,362,374,389,402,414,429,&
1192 &     453,468,480,490,500,510,520)
1193      brvlttbw=4     !A
1194    case(75,96,108,120,132,153,302,318,334,350,390,430,454)
1195      brvlttbw=5     !B
1196    case(5,14,23,29,63,76,84,97,109,121,133,141,154,163,&
1197 &     194,202,210,218,255,263,275,287,303,319,335,351,363,375,&
1198 &     391,403,415,431,439,455)
1199      brvlttbw=6     !C
1200    case(6,15,24,30,65,77,86,98,110,122,134,143,155,164,256,&
1201 &     264,276,288,304,320,336,352,364,376,392,404,416,432,440,456)
1202      brvlttbw=7     !I
1203    case(48,223,228,526,532)
1204      brvlttbw=8     !s
1205    end select
1206  case(75:142)       !Tetragonal
1207    select case(spgroupma)
1208    case(3,9,15,21,27,31,35,41,45:47,53:55,61:63,69:71,&
1209 &     77:79,83:85,89:91,97:99,105:107,113:115,121:123,129:131,&
1210 &     137:139,145:147,153:155,159:161,165:167,173:175,181:183,&
1211 &     189:191,197:199,205:207,213:215,221:223,229:231,235:237,&
1212 &     241:243,247:249,253:255,261:263,269:271,277:279,285:287,&
1213 &     293:295,301:303,309:311,317:319,323:325,329:331,335:337,&
1214 &     341:347,353:359,365:371,377:383,389:395,401:407,413:419,&
1215 &     425:431,437:443,449:455,461:467,473:479,485:491,497:503,&
1216 &     509:515,521:527,533:539,543:549,553:559,563:569)
1217      brvlttbw=9      !ShubIII
1218    case(4,10,16,22,28,32,36,42,48,56,64,72,80,86,92,100,&
1219 &     108,116,124,132,140,148,156,162,168,176,184,192,200,208,&
1220 &     216,224,232,238,244,250,256,264,272,280,288,296,304,312,&
1221 &     320,326,332,338,348,360,372,384,396,408,420,432,444,456,&
1222 &     468,480,492,504,516,528,540,550,560,570)
1223      brvlttbw=3     !c
1224    case(5,11,17,23,37,49,57,65,73,93,101,109,117,125,133,141,&
1225 &     149,169,177,185,193,201,209,217,225,257,265,273,281,289,297,&
1226 &     305,313,349,361,373,385,397,409,421,433,445,457,469,481,493,&
1227 &     505,517,529)
1228      brvlttbw=6     !C
1229    case(6,12,18,24,38,50,58,66,74,94,102,110,118,126,134,142,&
1230 &     150,170,178,186,194,202,210,218,226,258,266,274,282,290,298,&
1231 &     306,314,350,362,374,386,398,410,422,434,446,458,470,482,494,&
1232 &     506,518,530)
1233      brvlttbw=7     !I
1234    end select
1235  case(143:194)      !Hexagonal
1236    select case(spgroupma)
1237    case(15,19,23,27,31,35,39,43,47,51,55,59,63,67,71,&
1238 &     75:77,81:83,87:89,93:95,99:101,105:107,111,115,119,123,127,&
1239 &     131,135,139:141,145:147,151:153,157:159,163:165,169:171,&
1240 &     175:177,181:183,187:189,193:195,199:201,205:207,211:213,&
1241 &     217:219,223:225,229:231,235:241,245:251,255:261,265:271)
1242      brvlttbw=9      !ShubIII
1243    case(3,6,9,16,24,28,32,36,40,44,52,56,60,64,78,84,&
1244 &     90,96,112,116,120,124,128,132,136,142,148,154,160,166,172,&
1245 &     178,184,190,196,202,208,214,220,226,232,242,252,262,272)
1246      brvlttbw=6     !C
1247    case(12,20,48,68,72,102,108)
1248      brvlttbw=7     !I
1249    end select
1250  case(195:230)       !Cubic
1251    select case(spgroupma)
1252    case(16,20,24,28,32,35,39,42,46,50,54,58,61,65,69,&
1253 &     72,76,80,83,87,91,94:96,100:102,106:108,112:114,118:120,&
1254 &     124:126,130:132,136:138,142:144,147:149)
1255      brvlttbw=9      !ShubIII
1256    case(3,11,17,21,36,43,47,62,66,73,84,97,103,109,115)
1257      brvlttbw=7     !I
1258    case(6,25,29,51,55,77,88,121,127,133,139)
1259      brvlttbw=8     !s
1260    end select
1261  end select
1262 
1263  if(brvlttbw==9)shubnikov=3 !Shubnikov type III
1264  if(brvlttbw>=1 .and. brvlttbw<=8) shubnikov=4 !Shubnikov type IV
1265 
1266  genafm(:)=zero
1267  if(shubnikov==4)then
1268    if(brvlttbw==1)genafm(1)=half
1269    if(brvlttbw==2)genafm(2)=half
1270    if(brvlttbw==3)genafm(3)=half
1271    if(brvlttbw>=4 .and. brvlttbw<=8)then
1272      genafm(:)=half
1273      if(brvlttbw==4)genafm(1)=zero
1274      if(brvlttbw==5)genafm(2)=zero
1275      if(brvlttbw==6)genafm(3)=zero
1276    end if
1277  end if
1278 
1279 !DEBUG
1280 !write(std_out,*) 'gensymshub : end '
1281 !write(std_out,*) 'gensymshub, shubnikov =',shubnikov
1282 !ENDDEBUG
1283 
1284 end subroutine gensymshub

m_spgbuilder/gensymshub4 [ Functions ]

[ Top ] [ m_spgbuilder ] [ Functions ]

NAME

 gensymshub4

FUNCTION

 Assigns the Bravais magnetic translations for shubnikov type IV
 symmetry groups starting from the Fedorov space group
 and the translation, generator of the anti-ferromagnetic
 operations (as input). It will double nsym and tnons. In the end both
 symrel and tnons are completely determined.

INPUTS

 genafm(3) = translation, generator of the anti-ferromagnetic symmatry operations
 msym = maximum number of symmetry operations

OUTPUT

 symafm(msym)= (anti)ferromagnetic part of symmetry operations

SIDE EFFECTS

 nsym=number of symmetry operations, without magnetic operations at input,
  and with magnetic operations at output
 symrel(3,3,msym)=symmetry operations in real space in terms
  of primitive translations, without magnetic operations at input,
  and with magnetic operations at output
 tnons(3,msym)=nonsymmorphic translations for symmetry operations
  without magnetic operations at input,
  and with magnetic operations at output

PARENTS

      ingeo

CHILDREN

SOURCE

1322 subroutine gensymshub4(genafm,msym,nsym,symafm,symrel,tnons)
1323 
1324 
1325 !This section has been created automatically by the script Abilint (TD).
1326 !Do not modify the following lines by hand.
1327 #undef ABI_FUNC
1328 #define ABI_FUNC 'gensymshub4'
1329 !End of the abilint section
1330 
1331  implicit none
1332 
1333 !Arguments ------------------------------------
1334 !scalars
1335  integer,intent(in) :: msym
1336  integer,intent(inout) :: nsym
1337 !arrays
1338  integer,intent(inout) :: symafm(msym),symrel(3,3,msym)
1339  real(dp),intent(in) :: genafm(3)
1340  real(dp),intent(inout) :: tnons(3,msym)
1341 
1342 !Local variables ------------------------------
1343 !scalars
1344  integer :: ii
1345  character(len=500) :: message
1346 
1347 ! *************************************************************************
1348 
1349  if(msym<2*nsym)then
1350    write(message, '(3a)' )&
1351 &   'The number of symmetries in the Shubnikov type IV space group',ch10,&
1352 &   'is larger than the maximal allowed number of symmetries.'
1353    MSG_ERROR(message)
1354  end if
1355 
1356  do ii=1,nsym
1357    tnons(:,nsym+ii)=tnons(:,ii)+genafm(:)
1358    symrel(:,:,nsym+ii)=symrel(:,:,ii)
1359    symafm(ii)=1
1360    symafm(nsym+ii)=-1
1361  end do
1362  nsym=nsym*2
1363 
1364 end subroutine gensymshub4

m_spgbuilder/gensymspgr [ Functions ]

[ Top ] [ m_spgbuilder ] [ Functions ]

NAME

 gensymspgr

FUNCTION

 Give all the symmetry operations starting from the space group symbol.
 Suppose we are working in a conventional cell
 If brvltt 0 or positive, the pure translations of the
 Bravais lattice are included as generator of the group.
 In brvltt=-1, no pure translation is present, and the
 cell should be changed from conventional to primitive,
 outside of this routine.
 Treat also Shubnikov type III space groups.

INPUTS

 brvltt = input variable giving Bravais lattice
 msym = default number of symmetry operations
 shubnikov= magnetic type of the space group to be generated
 spgaxor = orientation of the cell axis (might be needed)
 spgorig = second choice of origin for certain groups
 (might be needed if nsym==0)
 spgroup = number of space group
 spgroupma= number of the magnetic space group

OUTPUT

 nsym = number of symmetry operations
 symafm(nsym)=(anti)ferromagnetic part of symmetry operations
 symrel(3,3,nsym)=symmetry operations in real space in terms
                  of primitive translations
 tnons(3,nsym)=nonsymmorphic translations for symmetry operations

SIDE EFFECTS

 brvltt = input variable giving Bravais lattice

PARENTS

      ingeo

CHILDREN

      chkgrp,print_symmetries,symsgcube,symsghexa,symsgmono,symsgortho
      symsgtetra

SOURCE

 92 subroutine gensymspgr(brvltt,msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma,symafm,symrel,tnons)
 93 
 94 
 95 !This section has been created automatically by the script Abilint (TD).
 96 !Do not modify the following lines by hand.
 97 #undef ABI_FUNC
 98 #define ABI_FUNC 'gensymspgr'
 99 !End of the abilint section
100 
101  implicit none
102 
103 !Arguments ------------------------------------
104 !scalars
105  integer,intent(in) :: msym,shubnikov,spgaxor,spgorig,spgroup,spgroupma
106  integer,intent(inout) :: brvltt
107  integer,intent(out) :: nsym
108 !arrays
109  integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i
110  real(dp),intent(inout) :: tnons(3,msym) !vz_i
111 
112 !Local variables ------------------------------
113 ! intsym,inttn = intermediate real to swap the columns of the symmetry matrix
114 ! bckbrvltt = backup to brvltt to compare the assigned and the input values
115 !real(dp) :: tsec(2)
116 !scalars
117  integer :: bckbrvltt,ii,inversion,jj,kk,ierr
118  integer :: intsym
119  real(dp) :: inttn
120  character(len=500) :: message
121 
122 ! *************************************************************************
123 
124 !List of the input parameters
125 !DEBUG
126 !write(std_out,*)' gensymspgr : enter with:'
127 !write(std_out,*)' spgroup = ',spgroup
128 !write(std_out,*)' spgaxor = ',spgaxor
129 !write(std_out,*)' spgorig = ',spgorig
130 !write(std_out,*)' brvltt  = ',brvltt
131 !ENDDEBUG
132 
133 !Assume that the value of spgroupma is consistent with the one of spgroup
134 !(this has been checked earlier)
135 
136 !Tests for consistency first the space group number and then the orientation
137 !Checks the space group number
138  if (.not.(spgroup>0 .and. spgroup<231) ) then
139    write(message, '(a,i0,a,a,a,a)' )&
140 &   'spgroup must be between 1 to 230, but is ',spgroup,ch10,&
141 &   'This is not allowed.  ',ch10,&
142 &   'Action: modify spgroup in the input file.'
143    MSG_ERROR(message)
144  end if
145 
146 !Checks the orientation
147  if (.not.(spgaxor>0 .and. spgaxor<10)) then
148    write(message, '(a,i12,a,a,a,a)' )&
149 &   'spgaxor must be from 1 to 9, but is',spgaxor,ch10,&
150 &   'This is not allowed.  ',ch10,&
151 &   'Action: modify spgaxor in the input file.'
152    MSG_ERROR(message)
153  end if
154 
155 !Checks the consistency between the origin and space group
156  if (spgorig==1 .or. spgorig==2) then
157  else
158    write(message, '(a,i4,a,a,a,a,a,a)' )&
159 &   'spgorig is',spgorig,ch10,&
160 &   'while it should be 0 or 1',ch10,&
161 &   'This is not allowed.  ',ch10,&
162 &   'Action: modify spgorig in the input file.'
163    MSG_ERROR(message)
164  end if
165 
166  if (spgorig>1) then
167    select case (spgroup)
168    case (48,50,59,68,70,85,86,88,125,126,129,130,133,134,137,138,&
169 &     141,142,201,203,222,224,227,228)
170    case default
171      write(message, '(a,a,a,a,a)' )&
172 &     'spgroup does not accept several origin choices',ch10,&
173 &     'This is not allowed.  ',ch10,&
174 &     'Action: modify spgorig in the input file.'
175      MSG_ERROR(message)
176    end select
177  end if
178 
179 !Checks for consistency between the orientation and space group
180  if (spgaxor>1) then
181    select case (spgroup)
182    case (3:74,146,148,155,160,161,166,167)
183    case default
184      write(message, '(a,a,a,a,a)' )&
185 &     'spgroup does not accept several orientations',ch10,&
186 &     'This is not allowed.  ',ch10,&
187 &     'Action: modify spgaxor or spgroup in the input file.'
188      MSG_ERROR(message)
189    end select
190  end if
191 
192  if (brvltt<-1 .or. brvltt>7)then
193    write(message, '(a,i4,a,a,a,a,a,a)' )&
194 &   'The input brvltt was ',brvltt,ch10,&
195 &   'and it should be an integer from -1 to 7',ch10,&
196 &   'This is not allowed.  ',ch10,&
197 &   'Action: modify brvltt  in the input file.'
198    MSG_ERROR(message)
199  end if
200 
201 !Assign nsym for each group according first to the order of the group
202 !Note that this value might be modified later:
203 !first because of the product with the inversion,
204 !second because of the centering operations
205  select case (spgroup)
206  case (1,2)
207    nsym=1
208  case (3:9)
209    nsym=2
210  case (143:148)
211    nsym=3
212  case (10:42,44:47,49,51:58,60:67,69,71:84,87)
213    nsym=4
214  case (149:176)
215    nsym=6
216  case (48,50,59,68,70,85,86,88:121,123,124,127,128,131,132,135,136,139,140)
217    nsym=8
218  case (177:200,202,204:206)
219    nsym=12
220  case (43,122,125,126,129,130,133,134,137,138)
221    nsym=16
222  case (201,207:219,221,223,225,226,229,230)
223    nsym=24
224  case (141,142)
225    nsym=32
226  case (203,220,222,224)
227    nsym=48
228  case (227,228)
229    nsym=192
230  end select
231 
232 !DEBUG
233 !write(std_out,*)'gensymspgr :  assigns nsym = ',nsym
234 !ENDDEBUG
235 
236 !Makes a backup to the brvltt for further comparison with the assigned value
237  bckbrvltt=brvltt
238 !Default brvltt
239  brvltt=1
240 
241 !call timab(47,1,tsec)
242 
243 !Assigns the first part of the symmetry operations:
244 !Rotation axis and mirror planes with or without translations,
245 !and sometimes also inversion operations.
246  select case (spgroup)
247  case (1:2)
248    symrel(:,:,1)=0
249    symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1
250    tnons(:,1)=zero ; symafm(1)=1
251  case (3:15)
252    call symsgmono(brvltt,msym,nsym,shubnikov,spgaxor,spgorig,spgroup,&
253 &   spgroupma,symafm,symrel,tnons)
254  case (16:74)
255    call symsgortho(msym,nsym,shubnikov,spgaxor,spgorig,spgroup,&
256 &   spgroupma,symafm,symrel,tnons)
257  case (75:142)
258    call symsgtetra(msym,nsym,shubnikov,spgaxor,spgorig,spgroup,&
259 &   spgroupma,symafm,symrel,tnons)
260  case (143:194)
261    call symsghexa(brvltt,msym,nsym,shubnikov,spgaxor,spgorig,spgroup,&
262 &   spgroupma,symafm,symrel,tnons)
263  case (195:230)
264    call symsgcube(msym,nsym,shubnikov,spgaxor,spgorig,spgroup,&
265 &   spgroupma,symafm,symrel,tnons)
266  end select
267 
268 !call timab(47,2,tsec)
269 
270 !Assign the inversion center (if necessary).
271 !Note that for monoclinic space groups, the inversion was already
272 !assigned in symsgmono.f. Some other inversions have also been assigned in the
273 !corresponding system routine
274  inversion=0
275  select case (spgroup)
276  case (2,47,49,51:58,60:67,69,71:74,83,84,87,123,124,127,128,131,132,&
277 &   135,136,139,140,147,148,162:167,175,176,191:194,200,202,204:206,&
278 &   221,223,225,226,229,230)
279    inversion=1
280 !    Treat the magnetic part
281    if(shubnikov==3)then
282      select case (spgroup)
283      case(2) ! Triclinic
284        inversion=-1
285      case(47,49,51:58,60:67,69,71:74) ! Orthorhombic
286        select case (spgroupma)
287        case(251,253,259,261,267,268,271,279,280,283,291,292,293,297,307,&
288 &         308,309,313,323,324,325,329,339,340,341,345,355,356,359,367,368,371,&
289 &         379,380,381,385,395,396,399,407,408,411,419,420,421,425,435,437,443,&
290 &         444,445,449,459,460,461,465,471,472,473,477,483,484,487,493,494,497,&
291 &         503,504,507,513,514,517,523,525,529,531,535,537,541,542,545,550,552,556,557,560)
292          inversion=-1
293        end select
294      case(83,84,87,123,124,127,128,131,132,135,136,139,140) ! Tetragonal
295        select case (spgroupma)
296        case(46,47,54,55,62,63,70,71,78,79,84,85,341,344,346,347,353,356,358,359,365,&
297 &         368,370,371,377,380,382,383,389,392,394,395,401,404,406,407,&
298 &         413,416,418,419,425,428,430,431,437,440,442,443,449,452,454,455,&
299 &         461,464,466,467,473,476,478,479,485,488,490,491,497,500,502,503,&
300 &         509,512,514,515,521,524,526,527,533,536,538,539,543,546,548,549)
301          inversion=-1
302        end select
303      case(147,148,162:167,175,176,191:194) ! Hexagonal or rhombohedral
304        select case (spgroupma)
305        case(15,19,75,76,81,82,87,88,93,94,99,100,105,106,139,140,145,146,&
306 &         235,236,237,241,245,246,247,251,255,256,257,261,265,266,267,271)
307          inversion=-1
308        end select
309      case(200,202,204:206,221,223,225,226,229,230) ! Cubic
310        select case (spgroupma)
311        case(16,20,24,28,32,35,39,94,96,100,102,106,108,112,114,118,120,124,126,&
312 &         130,132,136,138,142,144,147,149)
313          inversion=-1
314        end select
315      end select
316    end if
317  end select
318 
319 !DEBUG
320 !write(std_out,*)' gensymspgr : before inversion'
321 !write(std_out,*)' Describe the different symmetry operations (index,symrel,tnons,symafm)'
322 !do ii=1,nsym
323 !write(std_out,'(i3,2x,9i3,3es12.2,i3)')ii,symrel(:,:,ii),tnons(:,ii),symafm(ii)
324 !end do
325 !ENDDEBUG
326 
327  if(inversion/=0)then
328    do ii=1,nsym        ! visit all the symmetries assigned before
329      do jj=1,3        ! visit the 3x3 matrix corresponding to the symmetry i
330        tnons(jj,nsym+ii)=-tnons(jj,ii)
331        do kk=1,3
332          symrel(jj,kk,nsym+ii)=-symrel(jj,kk,ii)
333        end do
334      end do
335      symafm(nsym+ii)=inversion*symafm(ii)
336    end do
337    nsym=nsym*2
338  end if
339 
340 !DEBUG
341 !write(std_out,*)' gensymspgr : after inversion'
342 !write(std_out,*)' Describe the different symmetry operations (index,symrel,tnons,symafm)'
343 !do ii=1,nsym
344 !write(std_out,'(i3,2x,9i3,3es12.2,i3)')ii,symrel(:,:,ii),tnons(:,ii),symafm(ii)
345 !end do
346 !ENDDEBUG
347 
348 !Assign the Bravais lattice to each space group to which it has not yet
349 !been assigned
350  select case (spgroup)
351  case (38:41)
352    brvltt=5                ! A
353  case (20,21,35:37,63:68)
354    brvltt=4                ! C
355  case (22,42,43,69,70,196,202,203,209,210,216,219,225:228)
356    brvltt=3                ! F
357  case (23,24,44:46,71:74,79,80,82,87,88,97,98,107:110,119:122,&
358 &   139:142,197,199,204,206,211,214,217,220,229,230)
359    brvltt=2                ! I
360  case (146,148,155,160,161,166,167)
361    if (spgaxor==1) then
362      brvltt=7
363    end if
364  end select
365 
366  if (bckbrvltt/=0 .and. bckbrvltt/=-1) then
367    if (bckbrvltt/=brvltt) then
368      write(message, '(a,i8,a,a,a,i8,a,a)' )&
369 &     'The assigned brvltt ',brvltt,' is not equal',ch10,&
370 &     'to the input value ',bckbrvltt,ch10,&
371 &     'Assume experienced user. Execution will continue.'
372      MSG_WARNING(message)
373    end if
374  end if
375 
376 !if(bckbrvltt>=0)then
377 !Complete the set of primitive symmetries by translations
378 !associated with brvltt.
379  select case (brvltt)
380 !  Bravais lattice type : A ! translation associated: b/2+c/2
381  case (5)
382    do ii=1,nsym
383      tnons(1,nsym+ii)=tnons(1,ii)
384      tnons(2,nsym+ii)=tnons(2,ii)+0.5
385      tnons(3,nsym+ii)=tnons(3,ii)+0.5
386      symrel(:,:,nsym+ii)=symrel(:,:,ii)
387      symafm(nsym+ii)=symafm(ii)
388    end do
389    nsym=nsym*2
390 
391 !    Bravais lattice type : B ! translation associated: a/2+c/2
392  case (6)
393    do ii=1,nsym
394      tnons(1,nsym+ii)=tnons(1,ii)+0.5
395      tnons(2,nsym+ii)=tnons(2,ii)
396      tnons(3,nsym+ii)=tnons(3,ii)+0.5
397      symrel(:,:,nsym+ii)=symrel(:,:,ii)
398      symafm(nsym+ii)=symafm(ii)
399    end do
400    nsym=nsym*2
401 
402 !    Bravais lattice type : C ! translation associated: a/2+b/2
403  case (4)
404    do ii=1,nsym
405      tnons(1,nsym+ii)=tnons(1,ii)+0.5
406      tnons(2,nsym+ii)=tnons(2,ii)+0.5
407      tnons(3,nsym+ii)=tnons(3,ii)
408      symrel(:,:,nsym+ii)=symrel(:,:,ii)
409      symafm(nsym+ii)=symafm(ii)
410    end do
411    nsym=nsym*2
412 
413 !    Bravais lattice type : F ! translations associated: a/2+b/2,b/2+c/2,b/2+c/2
414  case (3)
415 !    For space groups containing d elements, all the symmetry operations
416 !    have already been obtained
417    if(spgroup/=43 .and. spgroup/=203 .and. spgroup/=227 .and. &
418 &   spgroup/=228)then
419      do ii=1,nsym
420 !        First translation: a/2+b/2
421        tnons(1,nsym+ii)=tnons(1,ii)+0.5
422        tnons(2,nsym+ii)=tnons(2,ii)+0.5
423        tnons(3,nsym+ii)=tnons(3,ii)
424        symrel(:,:,nsym+ii)=symrel(:,:,ii)
425        symafm(nsym+ii)=symafm(ii)
426      end do
427 !      Second translation: b/2+c/2
428      do ii=1,nsym
429        tnons(1,2*nsym+ii)=tnons(1,ii)
430        tnons(2,2*nsym+ii)=tnons(2,ii)+0.5
431        tnons(3,2*nsym+ii)=tnons(3,ii)+0.5
432        symrel(:,:,2*nsym+ii)=symrel(:,:,ii)
433        symafm(2*nsym+ii)=symafm(ii)
434      end do
435 !      Third translation: a/2+c/2
436      do ii=1,nsym
437        tnons(1,3*nsym+ii)=tnons(1,ii)+0.5
438        tnons(2,3*nsym+ii)=tnons(2,ii)
439        tnons(3,3*nsym+ii)=tnons(3,ii)+0.5
440        symrel(:,:,3*nsym+ii)=symrel(:,:,ii)
441        symafm(3*nsym+ii)=symafm(ii)
442      end do
443      nsym=nsym*4
444    end if
445 
446 !    Bravais lattice type: I ! translation associated: a/2+b/2+c/2
447  case (2)
448 !    For space groups containing d elements, all the symmetry operations
449 !    have already been obtained
450    if(spgroup/=122 .and. spgroup/=141 .and. spgroup/=142 .and. &
451 &   spgroup/=220 )then
452      do ii=1,nsym        ! visit all the symmetries assigned before
453        tnons(:,nsym+ii)=tnons(:,ii)+0.5
454        symrel(:,:,nsym+ii)=symrel(:,:,ii)
455        symafm(nsym+ii)=symafm(ii)
456      end do
457      nsym=nsym*2
458    end if
459 
460 !    Bravais lattice type: R
461 !    translations for hexagonal axes ONLY: (2/3,1/3,1/3) & (1/3,2/3,2/3)
462 !    first translation (2/3,1/3,1/3)
463  case (7)
464    do ii=1,nsym
465      tnons(1,nsym+ii)=tnons(1,ii)+two_thirds
466      tnons(2,nsym+ii)=tnons(2,ii)+third
467      tnons(3,nsym+ii)=tnons(3,ii)+third
468      symrel(:,:,nsym+ii)=symrel(:,:,ii)
469      symafm(nsym+ii)=symafm(ii)
470    end do
471 !    Second translation (1/3,2/3,2/3)
472    do ii=1,nsym
473      tnons(1,2*nsym+ii)=tnons(1,ii)+third
474      tnons(2,2*nsym+ii)=tnons(2,ii)+two_thirds
475      tnons(3,2*nsym+ii)=tnons(3,ii)+two_thirds
476      symrel(:,:,2*nsym+ii)=symrel(:,:,ii)
477      symafm(2*nsym+ii)=symafm(ii)
478    end do
479    nsym=nsym*3
480 
481  end select
482 
483 !end if
484 
485 !Translate tnons in the ]-0.5,0.5] interval
486  tnons(:,1:nsym)=tnons(:,1:nsym)-nint(tnons(:,1:nsym)-1.0d-8)
487 
488 !Orientations for the orthorhombic space groups
489 !WARNING : XG 000620 : I am not sure that this coding is correct !!
490  if (spgroup>15 .and. spgroup <75) then
491    select case (spgaxor)
492    case (1)             ! abc
493      write(std_out,*)' the choosen orientation corresponds to: abc; the proper one'
494    case (2)             ! cab
495      do ii=1,nsym
496        intsym=symrel(1,1,ii)
497        symrel(1,1,ii)=symrel(2,2,ii)
498        symrel(2,2,ii)=intsym
499        inttn=tnons(1,ii)
500        tnons(1,ii)=tnons(2,ii)
501        tnons(2,ii)=inttn
502      end do
503      do ii=1,nsym
504        intsym=symrel(1,1,ii)
505        symrel(1,1,ii)=symrel(3,3,ii)
506        symrel(3,3,ii)=intsym
507        inttn=tnons(1,ii)
508        tnons(1,ii)=tnons(3,ii)
509        tnons(3,ii)=inttn
510      end do
511      write(std_out,*)' the choosen orientation corresponds to:  cab'
512    case (3)             ! bca
513      do ii=1,nsym
514        intsym=symrel(1,1,ii)
515        symrel(1,1,ii)=symrel(2,2,ii)
516        symrel(2,2,ii)=intsym
517        inttn=tnons(1,ii)
518        tnons(1,ii)=tnons(2,ii)
519        tnons(2,ii)=inttn
520      end do
521      do ii=1,nsym
522        intsym=symrel(2,2,ii)
523        symrel(2,2,ii)=symrel(3,3,ii)
524        symrel(3,3,ii)=intsym
525        inttn=tnons(2,ii)
526        tnons(2,ii)=tnons(3,ii)
527        tnons(3,ii)=inttn
528      end do
529      write(std_out,*)' the choosen orientation corresponds to:  bca'
530    case (4)             ! acb
531      do ii=1,nsym
532        intsym=symrel(2,2,ii)
533        symrel(2,2,ii)=symrel(3,3,ii)
534        symrel(3,3,ii)=intsym
535        inttn=tnons(1,ii)
536        tnons(2,ii)=tnons(3,ii)
537        tnons(3,ii)=inttn
538      end do
539      write(std_out,*)' the choosen orientation corresponds to:  acb'
540    case (5)             ! bac
541      do ii=1,nsym
542        intsym=symrel(1,1,ii)
543        symrel(1,1,ii)=symrel(2,2,ii)
544        symrel(2,2,ii)=intsym
545        inttn=tnons(1,ii)
546        tnons(1,ii)=tnons(2,ii)
547        tnons(2,ii)=inttn
548      end do
549      write(std_out,*)' the choosen orientation corresponds to:  bac'
550    case (6)             ! cba
551      do ii=1,nsym
552        intsym=symrel(1,1,ii)
553        symrel(1,1,ii)=symrel(3,3,ii)
554        symrel(3,3,ii)=intsym
555        inttn=tnons(1,ii)
556        tnons(1,ii)=tnons(3,ii)
557        tnons(3,ii)=inttn
558      end do
559      write(std_out,*)' the choosen orientation corresponds to:  cba'
560    end select
561  end if
562 
563 !DEBUG
564 !write(std_out,*)' gensymspgr  : out of the Bravais lattice, nsym is',nsym
565 !ENDDEBUG
566 
567  call chkgrp(nsym,symafm,symrel,ierr)
568  if (ierr/=0) then
569    call print_symmetries(nsym,symrel,tnons,symafm)
570  end if
571 
572  ABI_CHECK(ierr==0,"Error in group closure")
573 
574 end subroutine gensymspgr