TABLE OF CONTENTS
ABINIT/m_spgbuilder [ Modules ]
NAME
m_spgbuilder
FUNCTION
Spacegroup builder.
COPYRIGHT
Copyright (C) 2008-2022 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 .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_spgbuilder 23 24 use defs_basis 25 use m_errors 26 use m_abicore 27 28 use m_symtk, only : chkgrp, print_symmetries 29 use m_symsg, only : symsgcube, symsghexa, symsgmono, symsgortho, symsgtetra 30 31 implicit none 32 33 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
SOURCE
581 subroutine gensymshub(genafm,spgroup,spgroupma,shubnikov) 582 583 !Arguments ------------------------------------ 584 !scalars 585 integer,intent(in) :: spgroup,spgroupma 586 integer,intent(out) :: shubnikov 587 !arrays 588 real(dp),intent(out) :: genafm(3) 589 590 !Local variables ------------------------------ 591 !scalars 592 integer :: brvlttbw=0,spgrmatch=1 593 character(len=500) :: message 594 595 ! ************************************************************************* 596 597 !List of the input parameters 598 !DEBUG 599 !write(std_out,*) 600 !write(std_out,*)' gensymshub : enter with:' 601 !write(std_out,*)' brvlttbw = ',brvlttbw 602 !write(std_out,*)' spgroup = ',spgroup 603 !write(std_out,*)' spgroupma = ',spgroupma 604 !ENDDEBUG 605 606 !Test for consistency the magnetic and non-magnetic space group 607 select case (spgroup) 608 case(1) 609 if (.not.(spgroupma>=3 .and. 3 >=spgroupma) ) spgrmatch=0 610 case(2) 611 if (.not.(spgroupma>=6 .and. 7 >=spgroupma) ) spgrmatch=0 612 case(3) 613 if (.not.(spgroupma>=3 .and. 6 >=spgroupma) ) spgrmatch=0 614 case(4) 615 if (.not.(spgroupma>=9 .and. 12 >=spgroupma) ) spgrmatch=0 616 case(5) 617 if (.not.(spgroupma>=15 .and. 17 >=spgroupma) ) spgrmatch=0 618 case(6) 619 if (.not.(spgroupma>=20 .and. 23 >=spgroupma) ) spgrmatch=0 620 case(7) 621 if (.not.(spgroupma>=26 .and. 31 >=spgroupma) ) spgrmatch=0 622 case(8) 623 if (.not.(spgroupma>=34 .and. 36 >=spgroupma) ) spgrmatch=0 624 case(9) 625 if (.not.(spgroupma>=39 .and. 41 >=spgroupma) ) spgrmatch=0 626 case(10) 627 if (.not.(spgroupma>=44 .and. 49 >=spgroupma) ) spgrmatch=0 628 case(11) 629 if (.not.(spgroupma>=52 .and. 57 >=spgroupma) ) spgrmatch=0 630 case(12) 631 if (.not.(spgroupma>=60 .and. 64 >=spgroupma) ) spgrmatch=0 632 case(13) 633 if (.not.(spgroupma>=67 .and. 74 >=spgroupma) ) spgrmatch=0 634 case(14) 635 if (.not.(spgroupma>=77 .and. 84 >=spgroupma) ) spgrmatch=0 636 case(15) 637 if (.not.(spgroupma>=87 .and. 91 >=spgroupma) ) spgrmatch=0 638 case(16) 639 if (.not.(spgroupma>= 3 .and. 6>=spgroupma) ) spgrmatch=0 640 case(17) 641 if (.not.(spgroupma>= 9 .and. 15 >=spgroupma) ) spgrmatch=0 642 case(18) 643 if (.not.(spgroupma>=18 .and. 24 >=spgroupma) ) spgrmatch=0 644 case(19) 645 if (.not.(spgroupma>=27 .and. 30 >=spgroupma) ) spgrmatch=0 646 case(20) 647 if (.not.(spgroupma>=33 .and. 37 >=spgroupma) ) spgrmatch=0 648 case(21) 649 if (.not.(spgroupma>=40 .and. 44 >=spgroupma) ) spgrmatch=0 650 case(22) 651 if (.not.(spgroupma>=47 .and. 48 >=spgroupma) ) spgrmatch=0 652 case(23) 653 if (.not.(spgroupma>=51 .and. 52 >=spgroupma) ) spgrmatch=0 654 case(24) 655 if (.not.(spgroupma>=55 .and. 56 >=spgroupma) ) spgrmatch=0 656 case(25) 657 if (.not.(spgroupma>=59 .and. 65 >=spgroupma) ) spgrmatch=0 658 case(26) 659 if (.not.(spgroupma>=68 .and. 77 >=spgroupma) ) spgrmatch=0 660 case(27) 661 if (.not.(spgroupma>=80 .and. 86 >=spgroupma) ) spgrmatch=0 662 case(28) 663 if (.not.(spgroupma>=89 .and. 98 >=spgroupma) ) spgrmatch=0 664 case(29) 665 if (.not.(spgroupma>=101 .and. 110 >=spgroupma) ) spgrmatch=0 666 case(30) 667 if (.not.(spgroupma>=113 .and. 122 >=spgroupma) ) spgrmatch=0 668 case(31) 669 if (.not.(spgroupma>=125 .and. 134 >=spgroupma) ) spgrmatch=0 670 case(32) 671 if (.not.(spgroupma>=137 .and. 143 >=spgroupma) ) spgrmatch=0 672 case(33) 673 if (.not.(spgroupma>=146 .and. 155 >=spgroupma) ) spgrmatch=0 674 case(34) 675 if (.not.(spgroupma>=158 .and. 164 >=spgroupma) ) spgrmatch=0 676 case(35) 677 if (.not.(spgroupma>=167 .and. 171>=spgroupma) ) spgrmatch=0 678 case(36) 679 if (.not.(spgroupma>=174 .and. 179>=spgroupma) ) spgrmatch=0 680 case(37) 681 if (.not.(spgroupma>=182 .and. 186 >=spgroupma) ) spgrmatch=0 682 case(38) 683 if (.not.(spgroupma>=189 .and. 194 >=spgroupma) ) spgrmatch=0 684 case(39) 685 if (.not.(spgroupma>=197 .and. 202 >=spgroupma) ) spgrmatch=0 686 case(40) 687 if (.not.(spgroupma>=205 .and. 210 >=spgroupma) ) spgrmatch=0 688 case(41) 689 if (.not.(spgroupma>=213 .and. 218 >=spgroupma) ) spgrmatch=0 690 case(42) 691 if (.not.(spgroupma>=221 .and. 223>=spgroupma) ) spgrmatch=0 692 case(43) 693 if (.not.(spgroupma>=226 .and. 228>=spgroupma) ) spgrmatch=0 694 case(44) 695 if (.not.(spgroupma>=231 .and. 234 >=spgroupma) ) spgrmatch=0 696 case(45) 697 if (.not.(spgroupma>=237 .and. 240 >=spgroupma) ) spgrmatch=0 698 case(46) 699 if (.not.(spgroupma>=243 .and. 248 >=spgroupma) ) spgrmatch=0 700 case(47) 701 if (.not.(spgroupma>=251 .and. 256 >=spgroupma) ) spgrmatch=0 702 case(48) 703 if (.not.(spgroupma>=259 .and. 264 >=spgroupma) ) spgrmatch=0 704 case(49) 705 if (.not.(spgroupma>=267 .and. 276 >=spgroupma) ) spgrmatch=0 706 case(50) 707 if (.not.(spgroupma>=279 .and. 288 >=spgroupma) ) spgrmatch=0 708 case(51) 709 if (.not.(spgroupma>=291 .and. 304 >=spgroupma) ) spgrmatch=0 710 case(52) 711 if (.not.(spgroupma>=307 .and. 320 >=spgroupma) ) spgrmatch=0 712 case(53) 713 if (.not.(spgroupma>=323 .and. 336 >=spgroupma) ) spgrmatch=0 714 case(54) 715 if (.not.(spgroupma>=339 .and. 352 >=spgroupma) ) spgrmatch=0 716 case(55) 717 if (.not.(spgroupma>=355 .and. 364 >=spgroupma) ) spgrmatch=0 718 case(56) 719 if (.not.(spgroupma>=367 .and. 376 >=spgroupma) ) spgrmatch=0 720 case(57) 721 if (.not.(spgroupma>=379 .and. 392 >=spgroupma) ) spgrmatch=0 722 case(58) 723 if (.not.(spgroupma>=395 .and. 404 >=spgroupma) ) spgrmatch=0 724 case(59) 725 if (.not.(spgroupma>=407 .and. 416 >=spgroupma) ) spgrmatch=0 726 case(60) 727 if (.not.(spgroupma>=419 .and. 432 >=spgroupma) ) spgrmatch=0 728 case(61) 729 if (.not.(spgroupma>=435 .and. 440 >=spgroupma) ) spgrmatch=0 730 case(62) 731 if (.not.(spgroupma>=443 .and. 456 >=spgroupma) ) spgrmatch=0 732 case(63) 733 if (.not.(spgroupma>=459 .and. 468 >=spgroupma) ) spgrmatch=0 734 case(64) 735 if (.not.(spgroupma>=471 .and. 480 >=spgroupma) ) spgrmatch=0 736 case(65) 737 if (.not.(spgroupma>=483 .and. 490 >=spgroupma) ) spgrmatch=0 738 case(66) 739 if (.not.(spgroupma>=493 .and. 500 >=spgroupma) ) spgrmatch=0 740 case(67) 741 if (.not.(spgroupma>=503 .and. 510 >=spgroupma) ) spgrmatch=0 742 case(68) 743 if (.not.(spgroupma>=513 .and. 520 >=spgroupma) ) spgrmatch=0 744 case(69) 745 if (.not.(spgroupma>=523 .and. 526 >=spgroupma) ) spgrmatch=0 746 case(70) 747 if (.not.(spgroupma>=529 .and. 532 >=spgroupma) ) spgrmatch=0 748 case(71) 749 if (.not.(spgroupma>=535 .and. 538 >=spgroupma) ) spgrmatch=0 750 case(72) 751 if (.not.(spgroupma>=541 .and. 547 >=spgroupma) ) spgrmatch=0 752 case(73) 753 if (.not.(spgroupma>=550 .and. 553 >=spgroupma) ) spgrmatch=0 754 case(74) 755 if (.not.(spgroupma>=556 .and. 562 >=spgroupma) ) spgrmatch=0 756 case(75) 757 if (.not.(spgroupma>= 3 .and. 6>=spgroupma) ) spgrmatch=0 758 case(76) 759 if (.not.(spgroupma>= 9 .and. 12 >=spgroupma) ) spgrmatch=0 760 case(77) 761 if (.not.(spgroupma>= 15 .and. 18>=spgroupma) ) spgrmatch=0 762 case(78) 763 if (.not.(spgroupma>= 21 .and. 24>=spgroupma) ) spgrmatch=0 764 case(79) 765 if (.not.(spgroupma>= 27 .and. 28 >=spgroupma) ) spgrmatch=0 766 case(80) 767 if (.not.(spgroupma>= 31 .and. 32 >=spgroupma) ) spgrmatch=0 768 case(81) 769 if (.not.(spgroupma>= 35 .and. 38 >=spgroupma) ) spgrmatch=0 770 case(82) 771 if (.not.(spgroupma>= 41 .and. 42 >=spgroupma) ) spgrmatch=0 772 case(83) 773 if (.not.(spgroupma>=45 .and. 50>=spgroupma) ) spgrmatch=0 774 case(84) 775 if (.not.(spgroupma>=53 .and. 58>=spgroupma) ) spgrmatch=0 776 case(85) 777 if (.not.(spgroupma>=61 .and. 66>=spgroupma) ) spgrmatch=0 778 case(86) 779 if (.not.(spgroupma>=69 .and. 74>=spgroupma) ) spgrmatch=0 780 case(87) 781 if (.not.(spgroupma>=77 .and. 80>=spgroupma) ) spgrmatch=0 782 case(88) 783 if (.not.(spgroupma>=83 .and. 86>=spgroupma) ) spgrmatch=0 784 case(89) 785 if (.not.(spgroupma>=89 .and. 94>=spgroupma) ) spgrmatch=0 786 case(90) 787 if (.not.(spgroupma>=97 .and. 102>=spgroupma) ) spgrmatch=0 788 case(91) 789 if (.not.(spgroupma>=105 .and. 110>=spgroupma) ) spgrmatch=0 790 case(92) 791 if (.not.(spgroupma>=113 .and. 118>=spgroupma) ) spgrmatch=0 792 case(93) 793 if (.not.(spgroupma>=121 .and. 126>=spgroupma) ) spgrmatch=0 794 case(94) 795 if (.not.(spgroupma>=129 .and. 134>=spgroupma) ) spgrmatch=0 796 case(95) 797 if (.not.(spgroupma>=137 .and. 142>=spgroupma) ) spgrmatch=0 798 case(96) 799 if (.not.(spgroupma>=145 .and. 150>=spgroupma) ) spgrmatch=0 800 case(97) 801 if (.not.(spgroupma>=153 .and. 156>=spgroupma) ) spgrmatch=0 802 case(98) 803 if (.not.(spgroupma>=159 .and. 162>=spgroupma) ) spgrmatch=0 804 case(99) 805 if (.not.(spgroupma>=165 .and. 170>=spgroupma) ) spgrmatch=0 806 case(100) 807 if (.not.(spgroupma>=173 .and. 178>=spgroupma) ) spgrmatch=0 808 case(101) 809 if (.not.(spgroupma>=181 .and. 186>=spgroupma) ) spgrmatch=0 810 case(102) 811 if (.not.(spgroupma>=189 .and. 194>=spgroupma) ) spgrmatch=0 812 case(103) 813 if (.not.(spgroupma>=197 .and. 202>=spgroupma) ) spgrmatch=0 814 case(104) 815 if (.not.(spgroupma>=205 .and. 210>=spgroupma) ) spgrmatch=0 816 case(105) 817 if (.not.(spgroupma>=213 .and. 218>=spgroupma) ) spgrmatch=0 818 case(106) 819 if (.not.(spgroupma>=221 .and. 226>=spgroupma) ) spgrmatch=0 820 case(107) 821 if (.not.(spgroupma>=229 .and. 232>=spgroupma) ) spgrmatch=0 822 case(108) 823 if (.not.(spgroupma>=235 .and. 238>=spgroupma) ) spgrmatch=0 824 case(109) 825 if (.not.(spgroupma>=241 .and. 244>=spgroupma) ) spgrmatch=0 826 case(110) 827 if (.not.(spgroupma>=247 .and. 250>=spgroupma) ) spgrmatch=0 828 case(111) 829 if (.not.(spgroupma>=253 .and. 258>=spgroupma) ) spgrmatch=0 830 case(112) 831 if (.not.(spgroupma>=261 .and. 266>=spgroupma) ) spgrmatch=0 832 case(113) 833 if (.not.(spgroupma>=269 .and. 274>=spgroupma) ) spgrmatch=0 834 case(114) 835 if (.not.(spgroupma>=277 .and. 282>=spgroupma) ) spgrmatch=0 836 case(115) 837 if (.not.(spgroupma>=285 .and. 290>=spgroupma) ) spgrmatch=0 838 case(116) 839 if (.not.(spgroupma>=293 .and. 298>=spgroupma) ) spgrmatch=0 840 case(117) 841 if (.not.(spgroupma>=301 .and. 306>=spgroupma) ) spgrmatch=0 842 case(118) 843 if (.not.(spgroupma>=309 .and. 314>=spgroupma) ) spgrmatch=0 844 case(119) 845 if (.not.(spgroupma>=317 .and. 320>=spgroupma) ) spgrmatch=0 846 case(120) 847 if (.not.(spgroupma>=323 .and. 326>=spgroupma) ) spgrmatch=0 848 case(121) 849 if (.not.(spgroupma>=329 .and. 332>=spgroupma) ) spgrmatch=0 850 case(122) 851 if (.not.(spgroupma>=335 .and. 338>=spgroupma) ) spgrmatch=0 852 case(123) 853 if (.not.(spgroupma>=341 .and. 350>=spgroupma) ) spgrmatch=0 854 case(124) 855 if (.not.(spgroupma>=353 .and. 362>=spgroupma) ) spgrmatch=0 856 case(125) 857 if (.not.(spgroupma>=365 .and. 374>=spgroupma) ) spgrmatch=0 858 case(126) 859 if (.not.(spgroupma>=377 .and. 386>=spgroupma) ) spgrmatch=0 860 case(127) 861 if (.not.(spgroupma>=389 .and. 398>=spgroupma) ) spgrmatch=0 862 case(128) 863 if (.not.(spgroupma>=401 .and. 410>=spgroupma) ) spgrmatch=0 864 case(129) 865 if (.not.(spgroupma>=413 .and. 422>=spgroupma) ) spgrmatch=0 866 case(130) 867 if (.not.(spgroupma>=425 .and. 434>=spgroupma) ) spgrmatch=0 868 case(131) 869 if (.not.(spgroupma>=437 .and. 446>=spgroupma) ) spgrmatch=0 870 case(132) 871 if (.not.(spgroupma>=449 .and. 458>=spgroupma) ) spgrmatch=0 872 case(133) 873 if (.not.(spgroupma>=461 .and. 470>=spgroupma) ) spgrmatch=0 874 case(134) 875 if (.not.(spgroupma>=473 .and. 482>=spgroupma) ) spgrmatch=0 876 case(135) 877 if (.not.(spgroupma>=485 .and. 494>=spgroupma) ) spgrmatch=0 878 case(136) 879 if (.not.(spgroupma>=497 .and. 506>=spgroupma) ) spgrmatch=0 880 case(137) 881 if (.not.(spgroupma>=509 .and. 518>=spgroupma) ) spgrmatch=0 882 case(138) 883 if (.not.(spgroupma>=521 .and. 530>=spgroupma) ) spgrmatch=0 884 case(139) 885 if (.not.(spgroupma>=533 .and. 540>=spgroupma) ) spgrmatch=0 886 case(140) 887 if (.not.(spgroupma>=543 .and. 550>=spgroupma) ) spgrmatch=0 888 case(141) 889 if (.not.(spgroupma>=553 .and. 560>=spgroupma) ) spgrmatch=0 890 case(142) 891 if (.not.(spgroupma>=563 .and. 570>=spgroupma) ) spgrmatch=0 892 case(143) 893 if (.not.(spgroupma>=3 .and. 3>=spgroupma) ) spgrmatch=0 894 case(144) 895 if (.not.(spgroupma>=6 .and. 6>=spgroupma) ) spgrmatch=0 896 case(145) 897 if (.not.(spgroupma>=9 .and. 9>=spgroupma) ) spgrmatch=0 898 case(146) 899 if (.not.(spgroupma>=12 .and. 12>=spgroupma) ) spgrmatch=0 900 case(147) 901 if (.not.(spgroupma>=15 .and. 16>=spgroupma) ) spgrmatch=0 902 case(148) 903 if (.not.(spgroupma>=19 .and. 20>=spgroupma) ) spgrmatch=0 904 case(149) 905 if (.not.(spgroupma>=23 .and. 24>=spgroupma) ) spgrmatch=0 906 case(150) 907 if (.not.(spgroupma>=27 .and. 28>=spgroupma) ) spgrmatch=0 908 case(151) 909 if (.not.(spgroupma>=31 .and. 32>=spgroupma) ) spgrmatch=0 910 case(152) 911 if (.not.(spgroupma>=35 .and. 36>=spgroupma) ) spgrmatch=0 912 case(153) 913 if (.not.(spgroupma>=39 .and. 40>=spgroupma) ) spgrmatch=0 914 case(154) 915 if (.not.(spgroupma>=43 .and. 44>=spgroupma) ) spgrmatch=0 916 case(155) 917 if (.not.(spgroupma>=47 .and. 48>=spgroupma) ) spgrmatch=0 918 case(156) 919 if (.not.(spgroupma>=51 .and. 52>=spgroupma) ) spgrmatch=0 920 case(157) 921 if (.not.(spgroupma>=55 .and. 56>=spgroupma) ) spgrmatch=0 922 case(158) 923 if (.not.(spgroupma>=59 .and. 60>=spgroupma) ) spgrmatch=0 924 case(159) 925 if (.not.(spgroupma>=63 .and. 64>=spgroupma) ) spgrmatch=0 926 case(160) 927 if (.not.(spgroupma>=67 .and. 68>=spgroupma) ) spgrmatch=0 928 case(161) 929 if (.not.(spgroupma>=71 .and. 72>=spgroupma) ) spgrmatch=0 930 case(162) 931 if (.not.(spgroupma>=75 .and. 78>=spgroupma) ) spgrmatch=0 932 case(163) 933 if (.not.(spgroupma>=81 .and. 84>=spgroupma) ) spgrmatch=0 934 case(164) 935 if (.not.(spgroupma>=87 .and. 90>=spgroupma) ) spgrmatch=0 936 case(165) 937 if (.not.(spgroupma>=93 .and. 96>=spgroupma) ) spgrmatch=0 938 case(166) 939 if (.not.(spgroupma>=99 .and. 102>=spgroupma) ) spgrmatch=0 940 case(167) 941 if (.not.(spgroupma>=105 .and. 108>=spgroupma) ) spgrmatch=0 942 case(168) 943 if (.not.(spgroupma>=111 .and. 112>=spgroupma) ) spgrmatch=0 944 case(169) 945 if (.not.(spgroupma>=115 .and. 116>=spgroupma) ) spgrmatch=0 946 case(170) 947 if (.not.(spgroupma>=119 .and. 120>=spgroupma) ) spgrmatch=0 948 case(171) 949 if (.not.(spgroupma>=123 .and. 124>=spgroupma) ) spgrmatch=0 950 case(172) 951 if (.not.(spgroupma>=127 .and. 128>=spgroupma) ) spgrmatch=0 952 case(173) 953 if (.not.(spgroupma>=131 .and. 132>=spgroupma) ) spgrmatch=0 954 case(174) 955 if (.not.(spgroupma>=135 .and. 136>=spgroupma) ) spgrmatch=0 956 case(175) 957 if (.not.(spgroupma>=139 .and. 142>=spgroupma) ) spgrmatch=0 958 case(176) 959 if (.not.(spgroupma>=145 .and. 148>=spgroupma) ) spgrmatch=0 960 case(177) 961 if (.not.(spgroupma>=151 .and. 154>=spgroupma) ) spgrmatch=0 962 case(178) 963 if (.not.(spgroupma>=157 .and. 160>=spgroupma) ) spgrmatch=0 964 case(179) 965 if (.not.(spgroupma>=163 .and. 166>=spgroupma) ) spgrmatch=0 966 case(180) 967 if (.not.(spgroupma>=169 .and. 172>=spgroupma) ) spgrmatch=0 968 case(181) 969 if (.not.(spgroupma>=175 .and. 178>=spgroupma) ) spgrmatch=0 970 case(182) 971 if (.not.(spgroupma>=181 .and. 184>=spgroupma) ) spgrmatch=0 972 case(183) 973 if (.not.(spgroupma>=187 .and. 190>=spgroupma) ) spgrmatch=0 974 case(184) 975 if (.not.(spgroupma>=193 .and. 196>=spgroupma) ) spgrmatch=0 976 case(185) 977 if (.not.(spgroupma>=199 .and. 202>=spgroupma) ) spgrmatch=0 978 case(186) 979 if (.not.(spgroupma>=205 .and. 208>=spgroupma) ) spgrmatch=0 980 case(187) 981 if (.not.(spgroupma>=211 .and. 214>=spgroupma) ) spgrmatch=0 982 case(188) 983 if (.not.(spgroupma>=217 .and. 220>=spgroupma) ) spgrmatch=0 984 case(189) 985 if (.not.(spgroupma>=223 .and. 226>=spgroupma) ) spgrmatch=0 986 case(190) 987 if (.not.(spgroupma>=229 .and. 232>=spgroupma) ) spgrmatch=0 988 case(191) 989 if (.not.(spgroupma>=235 .and. 242>=spgroupma) ) spgrmatch=0 990 case(192) 991 if (.not.(spgroupma>=245 .and. 252>=spgroupma) ) spgrmatch=0 992 case(193) 993 if (.not.(spgroupma>=255 .and. 262>=spgroupma) ) spgrmatch=0 994 case(194) 995 if (.not.(spgroupma>=265 .and. 272>=spgroupma) ) spgrmatch=0 996 case(195) 997 if (.not.(spgroupma>=3 .and. 3>=spgroupma) ) spgrmatch=0 998 case(196) 999 if (.not.(spgroupma>=6 .and. 6>=spgroupma) ) spgrmatch=0 1000 case(197) 1001 spgrmatch=0 1002 case(198) 1003 if (.not.(spgroupma>=11 .and. 11>=spgroupma) ) spgrmatch=0 1004 case(199) 1005 spgrmatch=0 1006 case(200) 1007 if (.not.(spgroupma>=16 .and. 17>=spgroupma) ) spgrmatch=0 1008 case(201) 1009 if (.not.(spgroupma>=20 .and. 21>=spgroupma) ) spgrmatch=0 1010 case(202) 1011 if (.not.(spgroupma>=24 .and. 25>=spgroupma) ) spgrmatch=0 1012 case(203) 1013 if (.not.(spgroupma>=28 .and. 29>=spgroupma) ) spgrmatch=0 1014 case(204) 1015 if (.not.(spgroupma>=32 .and. 32>=spgroupma) ) spgrmatch=0 1016 case(205) 1017 if (.not.(spgroupma>=35 .and. 36>=spgroupma) ) spgrmatch=0 1018 case(206) 1019 if (.not.(spgroupma>=39 .and. 39>=spgroupma) ) spgrmatch=0 1020 case(207) 1021 if (.not.(spgroupma>=42 .and. 43>=spgroupma) ) spgrmatch=0 1022 case(208) 1023 if (.not.(spgroupma>=46 .and. 47>=spgroupma) ) spgrmatch=0 1024 case(209) 1025 if (.not.(spgroupma>=50 .and. 51>=spgroupma) ) spgrmatch=0 1026 case(210) 1027 if (.not.(spgroupma>=54 .and. 55>=spgroupma) ) spgrmatch=0 1028 case(211) 1029 if (.not.(spgroupma>=58 .and. 58>=spgroupma) ) spgrmatch=0 1030 case(212) 1031 if (.not.(spgroupma>=61 .and. 62>=spgroupma) ) spgrmatch=0 1032 case(213) 1033 if (.not.(spgroupma>=65 .and. 66>=spgroupma) ) spgrmatch=0 1034 case(214) 1035 if (.not.(spgroupma>=69 .and. 69>=spgroupma) ) spgrmatch=0 1036 case(215) 1037 if (.not.(spgroupma>=72 .and. 73>=spgroupma) ) spgrmatch=0 1038 case(216) 1039 if (.not.(spgroupma>=76 .and. 77>=spgroupma) ) spgrmatch=0 1040 case(217) 1041 if (.not.(spgroupma>=80 .and. 80>=spgroupma) ) spgrmatch=0 1042 case(218) 1043 if (.not.(spgroupma>=83 .and. 84>=spgroupma) ) spgrmatch=0 1044 case(219) 1045 if (.not.(spgroupma>=87 .and. 88>=spgroupma) ) spgrmatch=0 1046 case(220) 1047 if (.not.(spgroupma>=91 .and. 91>=spgroupma) ) spgrmatch=0 1048 case(221) 1049 if (.not.(spgroupma>=94 .and. 97>=spgroupma) ) spgrmatch=0 1050 case(222) 1051 if (.not.(spgroupma>=100 .and. 103>=spgroupma) ) spgrmatch=0 1052 case(223) 1053 if (.not.(spgroupma>=106 .and. 109>=spgroupma) ) spgrmatch=0 1054 case(224) 1055 if (.not.(spgroupma>=112 .and. 115>=spgroupma) ) spgrmatch=0 1056 case(225) 1057 if (.not.(spgroupma>=118 .and. 121>=spgroupma) ) spgrmatch=0 1058 case(226) 1059 if (.not.(spgroupma>=124 .and. 127>=spgroupma) ) spgrmatch=0 1060 case(227) 1061 if (.not.(spgroupma>=130 .and. 133>=spgroupma) ) spgrmatch=0 1062 case(228) 1063 if (.not.(spgroupma>=136 .and. 139>=spgroupma) ) spgrmatch=0 1064 case(229) 1065 if (.not.(spgroupma>=142 .and. 144>=spgroupma) ) spgrmatch=0 1066 case(230) 1067 if (.not.(spgroupma>=147 .and. 149>=spgroupma) ) spgrmatch=0 1068 case default 1069 write(message, '(3a,i8,4a)' )& 1070 & 'The non-magnetic spacegroup is not specified ',ch10,& 1071 & 'while the magnetic space group is specified, spgroupma= ',spgroupma,ch10,& 1072 & 'This is not allowed. ',ch10,& 1073 & 'Action: specify spgroup in the input file.' 1074 ABI_ERROR(message) 1075 end select 1076 1077 if (spgrmatch==0) then 1078 write(message, '(a,i8,a,a,i8,4a)' )& 1079 & 'mismatch between the non-magnetic spacegroup ',spgroup,ch10,& 1080 & 'and the magnetic space group ',spgroupma,ch10,& 1081 & 'This is not allowed. ',ch10,& 1082 & 'Action: modify spgroup or spgroupma in the input file.' 1083 ABI_ERROR(message) 1084 end if 1085 1086 !DEBUG 1087 !write(std_out,*) ' gensymshub, after check the spgroup ... ' 1088 !ENDDEBUG 1089 1090 !Assign the magnetic Bravais lattice type from the magnetic space group number 1091 !As the magnetic space group number begins with 1 for EACH crystal system 1092 !we must first make our choice as a function of spgroup 1093 1094 !Convention : 1095 !brvlttbw = input variable giving Bravais black-and-white translation 1096 !(from 1 to 8 : Shubnikov type IV space group) 1097 !1,2,3 = 1/2 translation along a, b, or c respectively 1098 !4,5,6 = translation corresponding to A,B,C centering, respectively 1099 !7 = I centering 1100 !8 = (1/2 0 0) centering corresponding to normal F lattice 1101 !9 -> Shubnikov type III space group 1102 !Note that the use of the table 7.3 (p585) of Bradney and Cracknell 1103 !for the definition of the translation vectors 1104 !is extremely confusing, due to the strange choice of basis 1105 !vectors of table 3.1. 1106 !See table 7.4 (p588) for the spgroupma interpretation 1107 1108 select case(spgroup) 1109 case(1,2) !Triclinic 1110 select case(spgroupma) 1111 case(6) 1112 brvlttbw=9 !ShubIII 1113 case(3,7) 1114 brvlttbw=1 !Ps (note that it is not body centered, according to Table 1115 ! 7.3 of Bradley and Cracknell) 1116 end select 1117 case(3:15) !Monoclinic 1118 select case(spgroupma) 1119 case(3,9,15,20,26,34,39,44:46,52:54,60:62,67:69,77:79,87:89) 1120 brvlttbw=9 !ShubIII 1121 case(4,10,17,21,27,36,41,47,55,64,70,80,91) 1122 brvlttbw=1 !a 1123 case(5,11,22,29,48,56,71,81) 1124 brvlttbw=2 !b 1125 case(16,28,35,40,63,72,82,90) 1126 brvlttbw=3 !c 1127 case(31,73,83) 1128 brvlttbw=4 !A 1129 case(6,12,23,30,49,57,74,84) 1130 brvlttbw=6 !C 1131 end select 1132 case(16:74) !Orthorhombic 1133 select case(spgroupma) 1134 case(3,9,10,18,19,27,33,34,40,41,47,51,55,59,60,68:70,& 1135 & 80,81,89:91,101:103,113:115,125:127,137,138,146:148,158,159,& 1136 & 167,168,174:176,182,183,189:191,197:199,205:207,213:215,221,& 1137 & 222,226,227,231,232,237,238,243:245,251:253,259:261,267:271,& 1138 & 279:283,291:297,307:313,323:329,339:345,355:359,367:371,& 1139 & 379:385,395:399,407:411,419:425,435:437,443:449,459:465,& 1140 & 471:477,483:487,493:497,503:507,513:517,523:525,529:531,& 1141 & 535:537,541:545,550:552,556:560) 1142 brvlttbw=9 !ShubIII 1143 case(4,11,20,28,36,43,62,71,83,92,104,116,128,140,& 1144 & 149,160,170,178,185,192,200,208,216,234,240,247,254,262,272,& 1145 & 284,298,314,330,346,360,372,386,400,412,426,438,450,467,& 1146 & 479,489,499,509,519,547,562) 1147 brvlttbw=1 !a 1148 case(72,93,105,117,129,150,248,299,315,331,347,387,427,451) 1149 brvlttbw=2 !b 1150 case(12,21,35,42,52,56,61,73,82,94,106,118,130,139,& 1151 & 151,161,169,177,184,193,201,209,217,233,239,246,273,285,300,& 1152 & 316,332,348,361,373,388,401,413,428,452,466,478,488,498,& 1153 & 508,518,538,546,553,561) 1154 brvlttbw=3 !c 1155 case(13,22,37,44,64,74,85,95,107,119,131,142,152,162,& 1156 & 171,179,186,274,286,301,317,333,349,362,374,389,402,414,429,& 1157 & 453,468,480,490,500,510,520) 1158 brvlttbw=4 !A 1159 case(75,96,108,120,132,153,302,318,334,350,390,430,454) 1160 brvlttbw=5 !B 1161 case(5,14,23,29,63,76,84,97,109,121,133,141,154,163,& 1162 & 194,202,210,218,255,263,275,287,303,319,335,351,363,375,& 1163 & 391,403,415,431,439,455) 1164 brvlttbw=6 !C 1165 case(6,15,24,30,65,77,86,98,110,122,134,143,155,164,256,& 1166 & 264,276,288,304,320,336,352,364,376,392,404,416,432,440,456) 1167 brvlttbw=7 !I 1168 case(48,223,228,526,532) 1169 brvlttbw=8 !s 1170 end select 1171 case(75:142) !Tetragonal 1172 select case(spgroupma) 1173 case(3,9,15,21,27,31,35,41,45:47,53:55,61:63,69:71,& 1174 & 77:79,83:85,89:91,97:99,105:107,113:115,121:123,129:131,& 1175 & 137:139,145:147,153:155,159:161,165:167,173:175,181:183,& 1176 & 189:191,197:199,205:207,213:215,221:223,229:231,235:237,& 1177 & 241:243,247:249,253:255,261:263,269:271,277:279,285:287,& 1178 & 293:295,301:303,309:311,317:319,323:325,329:331,335:337,& 1179 & 341:347,353:359,365:371,377:383,389:395,401:407,413:419,& 1180 & 425:431,437:443,449:455,461:467,473:479,485:491,497:503,& 1181 & 509:515,521:527,533:539,543:549,553:559,563:569) 1182 brvlttbw=9 !ShubIII 1183 case(4,10,16,22,28,32,36,42,48,56,64,72,80,86,92,100,& 1184 & 108,116,124,132,140,148,156,162,168,176,184,192,200,208,& 1185 & 216,224,232,238,244,250,256,264,272,280,288,296,304,312,& 1186 & 320,326,332,338,348,360,372,384,396,408,420,432,444,456,& 1187 & 468,480,492,504,516,528,540,550,560,570) 1188 brvlttbw=3 !c 1189 case(5,11,17,23,37,49,57,65,73,93,101,109,117,125,133,141,& 1190 & 149,169,177,185,193,201,209,217,225,257,265,273,281,289,297,& 1191 & 305,313,349,361,373,385,397,409,421,433,445,457,469,481,493,& 1192 & 505,517,529) 1193 brvlttbw=6 !C 1194 case(6,12,18,24,38,50,58,66,74,94,102,110,118,126,134,142,& 1195 & 150,170,178,186,194,202,210,218,226,258,266,274,282,290,298,& 1196 & 306,314,350,362,374,386,398,410,422,434,446,458,470,482,494,& 1197 & 506,518,530) 1198 brvlttbw=7 !I 1199 end select 1200 case(143:194) !Hexagonal 1201 select case(spgroupma) 1202 case(15,19,23,27,31,35,39,43,47,51,55,59,63,67,71,& 1203 & 75:77,81:83,87:89,93:95,99:101,105:107,111,115,119,123,127,& 1204 & 131,135,139:141,145:147,151:153,157:159,163:165,169:171,& 1205 & 175:177,181:183,187:189,193:195,199:201,205:207,211:213,& 1206 & 217:219,223:225,229:231,235:241,245:251,255:261,265:271) 1207 brvlttbw=9 !ShubIII 1208 case(3,6,9,16,24,28,32,36,40,44,52,56,60,64,78,84,& 1209 & 90,96,112,116,120,124,128,132,136,142,148,154,160,166,172,& 1210 & 178,184,190,196,202,208,214,220,226,232,242,252,262,272) 1211 ! brvlttbw=6 !C XG230719 This is likely erroneous 1212 brvlttbw=3 !c 1213 case(12,20,48,68,72,102,108) 1214 brvlttbw=7 !I 1215 end select 1216 case(195:230) !Cubic 1217 select case(spgroupma) 1218 case(16,20,24,28,32,35,39,42,46,50,54,58,61,65,69,& 1219 & 72,76,80,83,87,91,94:96,100:102,106:108,112:114,118:120,& 1220 & 124:126,130:132,136:138,142:144,147:149) 1221 brvlttbw=9 !ShubIII 1222 case(3,11,17,21,36,43,47,62,66,73,84,97,103,109,115) 1223 brvlttbw=7 !I 1224 case(6,25,29,51,55,77,88,121,127,133,139) 1225 brvlttbw=8 !s 1226 end select 1227 end select 1228 1229 if(brvlttbw==9)shubnikov=3 !Shubnikov type III 1230 if(brvlttbw>=1 .and. brvlttbw<=8) shubnikov=4 !Shubnikov type IV 1231 1232 genafm(:)=zero 1233 if(shubnikov==4)then 1234 if(brvlttbw==1)genafm(1)=half 1235 if(brvlttbw==2)genafm(2)=half 1236 if(brvlttbw==3)genafm(3)=half 1237 if(brvlttbw>=4 .and. brvlttbw<=8)then 1238 genafm(:)=half 1239 if(brvlttbw==4)genafm(1)=zero 1240 if(brvlttbw==5)genafm(2)=zero 1241 if(brvlttbw==6)genafm(3)=zero 1242 end if 1243 end if 1244 1245 !DEBUG 1246 !write(std_out,*) 'gensymshub : end ' 1247 !write(std_out,*) 'gensymshub, shubnikov =',shubnikov 1248 !ENDDEBUG 1249 1250 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
SOURCE
1283 subroutine gensymshub4(genafm,msym,nsym,symafm,symrel,tnons) 1284 1285 !Arguments ------------------------------------ 1286 !scalars 1287 integer,intent(in) :: msym 1288 integer,intent(inout) :: nsym 1289 !arrays 1290 integer,intent(inout) :: symafm(msym),symrel(3,3,msym) 1291 real(dp),intent(in) :: genafm(3) 1292 real(dp),intent(inout) :: tnons(3,msym) 1293 1294 !Local variables ------------------------------ 1295 !scalars 1296 integer :: ii 1297 character(len=500) :: message 1298 1299 ! ************************************************************************* 1300 1301 if(msym<2*nsym)then 1302 write(message, '(3a)' )& 1303 & 'The number of symmetries in the Shubnikov type IV space group',ch10,& 1304 & 'is larger than the maximal allowed number of symmetries.' 1305 ABI_ERROR(message) 1306 end if 1307 1308 do ii=1,nsym 1309 tnons(:,nsym+ii)=tnons(:,ii)+genafm(:) 1310 symrel(:,:,nsym+ii)=symrel(:,:,ii) 1311 symafm(ii)=1 1312 symafm(nsym+ii)=-1 1313 end do 1314 nsym=nsym*2 1315 1316 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
SOURCE
80 subroutine gensymspgr(brvltt,msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma,symafm,symrel,tnons) 81 82 !Arguments ------------------------------------ 83 !scalars 84 integer,intent(in) :: msym,shubnikov,spgaxor,spgorig,spgroup,spgroupma 85 integer,intent(inout) :: brvltt 86 integer,intent(out) :: nsym 87 !arrays 88 integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i 89 real(dp),intent(inout) :: tnons(3,msym) !vz_i 90 91 !Local variables ------------------------------ 92 ! intsym,inttn = intermediate real to swap the columns of the symmetry matrix 93 ! bckbrvltt = backup to brvltt to compare the assigned and the input values 94 !real(dp) :: tsec(2) 95 !scalars 96 integer :: bckbrvltt,ii,inversion,jj,kk,ierr 97 integer :: intsym 98 real(dp) :: inttn 99 character(len=500) :: message 100 101 ! ************************************************************************* 102 103 !List of the input parameters 104 !DEBUG 105 !write(std_out,*)' gensymspgr : enter with:' 106 !write(std_out,*)' spgroup = ',spgroup 107 !write(std_out,*)' spgaxor = ',spgaxor 108 !write(std_out,*)' spgorig = ',spgorig 109 !write(std_out,*)' brvltt = ',brvltt 110 !ENDDEBUG 111 112 !Assume that the value of spgroupma is consistent with the one of spgroup 113 !(this has been checked earlier) 114 115 !Tests for consistency first the space group number and then the orientation 116 !Checks the space group number 117 if (.not.(spgroup>0 .and. spgroup<231) ) then 118 write(message, '(a,i0,a,a,a,a)' )& 119 & 'spgroup must be between 1 to 230, but is ',spgroup,ch10,& 120 & 'This is not allowed. ',ch10,& 121 & 'Action: modify spgroup in the input file.' 122 ABI_ERROR(message) 123 end if 124 125 !Checks the orientation 126 if (.not.(spgaxor>0 .and. spgaxor<10)) then 127 write(message, '(a,i12,a,a,a,a)' )& 128 & 'spgaxor must be from 1 to 9, but is',spgaxor,ch10,& 129 & 'This is not allowed. ',ch10,& 130 & 'Action: modify spgaxor in the input file.' 131 ABI_ERROR(message) 132 end if 133 134 !Checks the consistency between the origin and space group 135 if (spgorig==1 .or. spgorig==2) then 136 else 137 write(message, '(a,i4,a,a,a,a,a,a)' )& 138 & 'spgorig is',spgorig,ch10,& 139 & 'while it should be 0 or 1',ch10,& 140 & 'This is not allowed. ',ch10,& 141 & 'Action: modify spgorig in the input file.' 142 ABI_ERROR(message) 143 end if 144 145 if (spgorig>1) then 146 select case (spgroup) 147 case (48,50,59,68,70,85,86,88,125,126,129,130,133,134,137,138,& 148 & 141,142,201,203,222,224,227,228) 149 case default 150 write(message, '(a,a,a,a,a)' )& 151 & 'spgroup does not accept several origin choices',ch10,& 152 & 'This is not allowed. ',ch10,& 153 & 'Action: modify spgorig in the input file.' 154 ABI_ERROR(message) 155 end select 156 end if 157 158 !Checks for consistency between the orientation and space group 159 if (spgaxor>1) then 160 select case (spgroup) 161 case (3:74,146,148,155,160,161,166,167) 162 case default 163 write(message, '(a,a,a,a,a)' )& 164 & 'spgroup does not accept several orientations',ch10,& 165 & 'This is not allowed. ',ch10,& 166 & 'Action: modify spgaxor or spgroup in the input file.' 167 ABI_ERROR(message) 168 end select 169 end if 170 171 if (brvltt<-1 .or. brvltt>7)then 172 write(message, '(a,i4,a,a,a,a,a,a)' )& 173 & 'The input brvltt was ',brvltt,ch10,& 174 & 'and it should be an integer from -1 to 7',ch10,& 175 & 'This is not allowed. ',ch10,& 176 & 'Action: modify brvltt in the input file.' 177 ABI_ERROR(message) 178 end if 179 180 !Assign nsym for each group according first to the order of the group 181 !Note that this value might be modified later: 182 !first because of the product with the inversion, 183 !second because of the centering operations 184 select case (spgroup) 185 case (1,2) 186 nsym=1 187 case (3:9) 188 nsym=2 189 case (143:148) 190 nsym=3 191 case (10:42,44:47,49,51:58,60:67,69,71:84,87) 192 nsym=4 193 case (149:176) 194 nsym=6 195 case (48,50,59,68,70,85,86,88:121,123,124,127,128,131,132,135,136,139,140) 196 nsym=8 197 case (177:200,202,204:206) 198 nsym=12 199 case (43,122,125,126,129,130,133,134,137,138) 200 nsym=16 201 case (201,207:219,221,223,225,226,229,230) 202 nsym=24 203 case (141,142) 204 nsym=32 205 case (203,220,222,224) 206 nsym=48 207 case (227,228) 208 nsym=192 209 end select 210 211 !DEBUG 212 !write(std_out,*)'gensymspgr : assigns nsym = ',nsym 213 !ENDDEBUG 214 215 !Makes a backup to the brvltt for further comparison with the assigned value 216 bckbrvltt=brvltt 217 !Default brvltt 218 brvltt=1 219 220 !call timab(47,1,tsec) 221 222 !Assigns the first part of the symmetry operations: 223 !Rotation axis and mirror planes with or without translations, 224 !and sometimes also inversion operations. 225 select case (spgroup) 226 case (1:2) 227 symrel(:,:,1)=0 228 symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 229 tnons(:,1)=zero ; symafm(1)=1 230 case (3:15) 231 call symsgmono(brvltt,msym,nsym,shubnikov,spgaxor,spgorig,spgroup,& 232 & spgroupma,symafm,symrel,tnons) 233 case (16:74) 234 call symsgortho(msym,nsym,shubnikov,spgaxor,spgorig,spgroup,& 235 & spgroupma,symafm,symrel,tnons) 236 case (75:142) 237 call symsgtetra(msym,nsym,shubnikov,spgaxor,spgorig,spgroup,& 238 & spgroupma,symafm,symrel,tnons) 239 case (143:194) 240 call symsghexa(brvltt,msym,nsym,shubnikov,spgaxor,spgorig,spgroup,& 241 & spgroupma,symafm,symrel,tnons) 242 case (195:230) 243 call symsgcube(msym,nsym,shubnikov,spgaxor,spgorig,spgroup,& 244 & spgroupma,symafm,symrel,tnons) 245 end select 246 247 !call timab(47,2,tsec) 248 249 !Assign the inversion center (if necessary). 250 !Note that for monoclinic space groups, the inversion was already 251 !assigned in symsgmono.f. Some other inversions have also been assigned in the 252 !corresponding system routine 253 inversion=0 254 select case (spgroup) 255 case (2,47,49,51:58,60:67,69,71:74,83,84,87,123,124,127,128,131,132,& 256 & 135,136,139,140,147,148,162:167,175,176,191:194,200,202,204:206,& 257 & 221,223,225,226,229,230) 258 inversion=1 259 ! Treat the magnetic part 260 if(shubnikov==3)then 261 select case (spgroup) 262 case(2) ! Triclinic 263 inversion=-1 264 case(47,49,51:58,60:67,69,71:74) ! Orthorhombic 265 select case (spgroupma) 266 case(251,253,259,261,267,268,271,279,280,283,291,292,293,297,307,& 267 & 308,309,313,323,324,325,329,339,340,341,345,355,356,359,367,368,371,& 268 & 379,380,381,385,395,396,399,407,408,411,419,420,421,425,435,437,443,& 269 & 444,445,449,459,460,461,465,471,472,473,477,483,484,487,493,494,497,& 270 & 503,504,507,513,514,517,523,525,529,531,535,537,541,542,545,550,552,556,557,560) 271 inversion=-1 272 end select 273 case(83,84,87,123,124,127,128,131,132,135,136,139,140) ! Tetragonal 274 select case (spgroupma) 275 case(46,47,54,55,62,63,70,71,78,79,84,85,341,344,346,347,353,356,358,359,365,& 276 & 368,370,371,377,380,382,383,389,392,394,395,401,404,406,407,& 277 & 413,416,418,419,425,428,430,431,437,440,442,443,449,452,454,455,& 278 & 461,464,466,467,473,476,478,479,485,488,490,491,497,500,502,503,& 279 & 509,512,514,515,521,524,526,527,533,536,538,539,543,546,548,549) 280 inversion=-1 281 end select 282 case(147,148,162:167,175,176,191:194) ! Hexagonal or rhombohedral 283 select case (spgroupma) 284 case(15,19,75,76,81,82,87,88,93,94,99,100,105,106,139,140,145,146,& 285 & 235,236,237,241,245,246,247,251,255,256,257,261,265,266,267,271) 286 inversion=-1 287 end select 288 case(200,202,204:206,221,223,225,226,229,230) ! Cubic 289 select case (spgroupma) 290 case(16,20,24,28,32,35,39,94,96,100,102,106,108,112,114,118,120,124,126,& 291 & 130,132,136,138,142,144,147,149) 292 inversion=-1 293 end select 294 end select 295 end if 296 end select 297 298 !DEBUG 299 !write(std_out,*)' gensymspgr : before inversion' 300 !write(std_out,*)' Describe the different symmetry operations (index,symrel,tnons,symafm)' 301 !do ii=1,nsym 302 !write(std_out,'(i3,2x,9i3,3es12.2,i3)')ii,symrel(:,:,ii),tnons(:,ii),symafm(ii) 303 !end do 304 !ENDDEBUG 305 306 if(inversion/=0)then 307 do ii=1,nsym ! visit all the symmetries assigned before 308 do jj=1,3 ! visit the 3x3 matrix corresponding to the symmetry i 309 tnons(jj,nsym+ii)=-tnons(jj,ii) 310 do kk=1,3 311 symrel(jj,kk,nsym+ii)=-symrel(jj,kk,ii) 312 end do 313 end do 314 symafm(nsym+ii)=inversion*symafm(ii) 315 end do 316 nsym=nsym*2 317 end if 318 319 !DEBUG 320 !write(std_out,*)' gensymspgr : after 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 !Assign the Bravais lattice to each space group to which it has not yet 328 !been assigned 329 select case (spgroup) 330 case (38:41) 331 brvltt=5 ! A 332 case (20,21,35:37,63:68) 333 brvltt=4 ! C 334 case (22,42,43,69,70,196,202,203,209,210,216,219,225:228) 335 brvltt=3 ! F 336 case (23,24,44:46,71:74,79,80,82,87,88,97,98,107:110,119:122,& 337 & 139:142,197,199,204,206,211,214,217,220,229,230) 338 brvltt=2 ! I 339 case (146,148,155,160,161,166,167) 340 if (spgaxor==1) then 341 brvltt=7 342 end if 343 end select 344 345 if (bckbrvltt/=0 .and. bckbrvltt/=-1) then 346 if (bckbrvltt/=brvltt) then 347 write(message, '(a,i8,a,a,a,i8,a,a)' )& 348 & 'The assigned brvltt ',brvltt,' is not equal',ch10,& 349 & 'to the input value ',bckbrvltt,ch10,& 350 & 'Assume experienced user. Execution will continue.' 351 ABI_WARNING(message) 352 end if 353 end if 354 355 !if(bckbrvltt>=0)then 356 !Complete the set of primitive symmetries by translations 357 !associated with brvltt. 358 select case (brvltt) 359 ! Bravais lattice type : A ! translation associated: b/2+c/2 360 case (5) 361 do ii=1,nsym 362 tnons(1,nsym+ii)=tnons(1,ii) 363 tnons(2,nsym+ii)=tnons(2,ii)+0.5 364 tnons(3,nsym+ii)=tnons(3,ii)+0.5 365 symrel(:,:,nsym+ii)=symrel(:,:,ii) 366 symafm(nsym+ii)=symafm(ii) 367 end do 368 nsym=nsym*2 369 370 ! Bravais lattice type : B ! translation associated: a/2+c/2 371 case (6) 372 do ii=1,nsym 373 tnons(1,nsym+ii)=tnons(1,ii)+0.5 374 tnons(2,nsym+ii)=tnons(2,ii) 375 tnons(3,nsym+ii)=tnons(3,ii)+0.5 376 symrel(:,:,nsym+ii)=symrel(:,:,ii) 377 symafm(nsym+ii)=symafm(ii) 378 end do 379 nsym=nsym*2 380 381 ! Bravais lattice type : C ! translation associated: a/2+b/2 382 case (4) 383 do ii=1,nsym 384 tnons(1,nsym+ii)=tnons(1,ii)+0.5 385 tnons(2,nsym+ii)=tnons(2,ii)+0.5 386 tnons(3,nsym+ii)=tnons(3,ii) 387 symrel(:,:,nsym+ii)=symrel(:,:,ii) 388 symafm(nsym+ii)=symafm(ii) 389 end do 390 nsym=nsym*2 391 392 ! Bravais lattice type : F ! translations associated: a/2+b/2,b/2+c/2,b/2+c/2 393 case (3) 394 ! For space groups containing d elements, all the symmetry operations 395 ! have already been obtained 396 if(spgroup/=43 .and. spgroup/=203 .and. spgroup/=227 .and. & 397 & spgroup/=228)then 398 do ii=1,nsym 399 ! First translation: a/2+b/2 400 tnons(1,nsym+ii)=tnons(1,ii)+0.5 401 tnons(2,nsym+ii)=tnons(2,ii)+0.5 402 tnons(3,nsym+ii)=tnons(3,ii) 403 symrel(:,:,nsym+ii)=symrel(:,:,ii) 404 symafm(nsym+ii)=symafm(ii) 405 end do 406 ! Second translation: b/2+c/2 407 do ii=1,nsym 408 tnons(1,2*nsym+ii)=tnons(1,ii) 409 tnons(2,2*nsym+ii)=tnons(2,ii)+0.5 410 tnons(3,2*nsym+ii)=tnons(3,ii)+0.5 411 symrel(:,:,2*nsym+ii)=symrel(:,:,ii) 412 symafm(2*nsym+ii)=symafm(ii) 413 end do 414 ! Third translation: a/2+c/2 415 do ii=1,nsym 416 tnons(1,3*nsym+ii)=tnons(1,ii)+0.5 417 tnons(2,3*nsym+ii)=tnons(2,ii) 418 tnons(3,3*nsym+ii)=tnons(3,ii)+0.5 419 symrel(:,:,3*nsym+ii)=symrel(:,:,ii) 420 symafm(3*nsym+ii)=symafm(ii) 421 end do 422 nsym=nsym*4 423 end if 424 425 ! Bravais lattice type: I ! translation associated: a/2+b/2+c/2 426 case (2) 427 ! For space groups containing d elements, all the symmetry operations 428 ! have already been obtained 429 if(spgroup/=122 .and. spgroup/=141 .and. spgroup/=142 .and. & 430 & spgroup/=220 )then 431 do ii=1,nsym ! visit all the symmetries assigned before 432 tnons(:,nsym+ii)=tnons(:,ii)+0.5 433 symrel(:,:,nsym+ii)=symrel(:,:,ii) 434 symafm(nsym+ii)=symafm(ii) 435 end do 436 nsym=nsym*2 437 end if 438 439 ! Bravais lattice type: R 440 ! translations for hexagonal axes ONLY: (2/3,1/3,1/3) & (1/3,2/3,2/3) 441 ! first translation (2/3,1/3,1/3) 442 case (7) 443 do ii=1,nsym 444 tnons(1,nsym+ii)=tnons(1,ii)+two_thirds 445 tnons(2,nsym+ii)=tnons(2,ii)+third 446 tnons(3,nsym+ii)=tnons(3,ii)+third 447 symrel(:,:,nsym+ii)=symrel(:,:,ii) 448 symafm(nsym+ii)=symafm(ii) 449 end do 450 ! Second translation (1/3,2/3,2/3) 451 do ii=1,nsym 452 tnons(1,2*nsym+ii)=tnons(1,ii)+third 453 tnons(2,2*nsym+ii)=tnons(2,ii)+two_thirds 454 tnons(3,2*nsym+ii)=tnons(3,ii)+two_thirds 455 symrel(:,:,2*nsym+ii)=symrel(:,:,ii) 456 symafm(2*nsym+ii)=symafm(ii) 457 end do 458 nsym=nsym*3 459 460 end select 461 462 !end if 463 464 !Translate tnons in the ]-0.5,0.5] interval 465 tnons(:,1:nsym)=tnons(:,1:nsym)-nint(tnons(:,1:nsym)-1.0d-8) 466 467 !Orientations for the orthorhombic space groups 468 !WARNING : XG 000620 : I am not sure that this coding is correct !! 469 if (spgroup>15 .and. spgroup <75) then 470 select case (spgaxor) 471 case (1) ! abc 472 write(std_out,*)' the choosen orientation corresponds to: abc; the proper one' 473 case (2) ! cab 474 do ii=1,nsym 475 intsym=symrel(1,1,ii) 476 symrel(1,1,ii)=symrel(2,2,ii) 477 symrel(2,2,ii)=intsym 478 inttn=tnons(1,ii) 479 tnons(1,ii)=tnons(2,ii) 480 tnons(2,ii)=inttn 481 end do 482 do ii=1,nsym 483 intsym=symrel(1,1,ii) 484 symrel(1,1,ii)=symrel(3,3,ii) 485 symrel(3,3,ii)=intsym 486 inttn=tnons(1,ii) 487 tnons(1,ii)=tnons(3,ii) 488 tnons(3,ii)=inttn 489 end do 490 write(std_out,*)' the choosen orientation corresponds to: cab' 491 case (3) ! bca 492 do ii=1,nsym 493 intsym=symrel(1,1,ii) 494 symrel(1,1,ii)=symrel(2,2,ii) 495 symrel(2,2,ii)=intsym 496 inttn=tnons(1,ii) 497 tnons(1,ii)=tnons(2,ii) 498 tnons(2,ii)=inttn 499 end do 500 do ii=1,nsym 501 intsym=symrel(2,2,ii) 502 symrel(2,2,ii)=symrel(3,3,ii) 503 symrel(3,3,ii)=intsym 504 inttn=tnons(2,ii) 505 tnons(2,ii)=tnons(3,ii) 506 tnons(3,ii)=inttn 507 end do 508 write(std_out,*)' the choosen orientation corresponds to: bca' 509 case (4) ! acb 510 do ii=1,nsym 511 intsym=symrel(2,2,ii) 512 symrel(2,2,ii)=symrel(3,3,ii) 513 symrel(3,3,ii)=intsym 514 inttn=tnons(1,ii) 515 tnons(2,ii)=tnons(3,ii) 516 tnons(3,ii)=inttn 517 end do 518 write(std_out,*)' the choosen orientation corresponds to: acb' 519 case (5) ! bac 520 do ii=1,nsym 521 intsym=symrel(1,1,ii) 522 symrel(1,1,ii)=symrel(2,2,ii) 523 symrel(2,2,ii)=intsym 524 inttn=tnons(1,ii) 525 tnons(1,ii)=tnons(2,ii) 526 tnons(2,ii)=inttn 527 end do 528 write(std_out,*)' the choosen orientation corresponds to: bac' 529 case (6) ! cba 530 do ii=1,nsym 531 intsym=symrel(1,1,ii) 532 symrel(1,1,ii)=symrel(3,3,ii) 533 symrel(3,3,ii)=intsym 534 inttn=tnons(1,ii) 535 tnons(1,ii)=tnons(3,ii) 536 tnons(3,ii)=inttn 537 end do 538 write(std_out,*)' the choosen orientation corresponds to: cba' 539 end select 540 end if 541 542 !DEBUG 543 !write(std_out,*)' gensymspgr : out of the Bravais lattice, nsym is',nsym 544 !ENDDEBUG 545 546 call chkgrp(nsym,symafm,symrel,ierr) 547 if (ierr/=0) then 548 call print_symmetries(nsym,symrel,tnons,symafm) 549 end if 550 551 ABI_CHECK(ierr==0,"Error in group closure") 552 553 end subroutine gensymspgr