TABLE OF CONTENTS
ABINIT/symlist_others [ Functions ]
NAME
symlist_others
FUNCTION
Determine the space group from the number and type of symmetry operations Non primitive, non BCC, non FCC case : rhombohedral or face centered
COPYRIGHT
Copyright (C) 2000-2018 ABINIT group (RC) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
brvltt=Bravais lattice type nsym=actual number of symmetries n_axes(31)=array containing the number of all the possible symmetry operations
OUTPUT
spgroup=space group number ; returns 0 if not found
NOTES
The list of symmetry operations is for the conventional cell
TODO
For the time being there are several groups where uncertainties still exist This will be solved in the very next ABINIT version
PARENTS
symspgr
CHILDREN
SOURCE
40 #if defined HAVE_CONFIG_H 41 #include "config.h" 42 #endif 43 44 #include "abi_common.h" 45 46 47 subroutine symlist_others(brvltt,nsym,n_axes,spgroup) 48 49 use defs_basis 50 use m_profiling_abi 51 52 !This section has been created automatically by the script Abilint (TD). 53 !Do not modify the following lines by hand. 54 #undef ABI_FUNC 55 #define ABI_FUNC 'symlist_others' 56 !End of the abilint section 57 58 implicit none 59 60 !Arguments ------------------------------------ 61 !scalars 62 integer,intent(in) :: brvltt,nsym 63 integer,intent(out) :: spgroup 64 !arrays 65 integer,intent(in) :: n_axes(31) 66 67 !Local variables------------------------------- 68 !character(len=500) :: message 69 !arrays 70 integer :: n_axest(31) 71 72 !************************************************************************** 73 74 !DEBUG 75 !write(std_out,*) ' symlist_others : enter ' 76 !write(std_out,*) ' nsym = ', nsym 77 !write(std_out,*) ' brvltt = ',brvltt 78 !write(std_out,'(a,10i3)' ) ' n_axes(1:10) =',n_axes(1:10) 79 !write(std_out,'(a,10i3)' ) ' n_axes(11:20)=',n_axes(11:20) 80 !write(std_out,'(a,11i3)' ) ' n_axes(21:31)=',n_axes(21:31) 81 !ENDDEBUG 82 83 spgroup=0 84 85 if(brvltt==4 .or. brvltt==5 .or. brvltt==6)then ! A, B, C face centered 86 87 select case(nsym) 88 89 case(4) 90 91 n_axest=(/0,0,0,0,0,0,1,1,1,0, 0,0,0,0,0,0,0,0,0,1, 0,0,0,0,0,0,0,0,0,0,0/) 92 if(sum((n_axes-n_axest)**2)==0) spgroup=5 93 n_axest=(/0,0,0,0,0,0,1,1,0,0, 0,0,0,0,1,1,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0/) 94 if(sum((n_axes-n_axest)**2)==0) spgroup=8 95 n_axest=(/0,0,0,0,0,0,1,1,0,0, 0,0,0,0,0,1,0,1,0,0, 0,0,0,0,0,0,0,0,0,0,0/) 96 if(sum((n_axes-n_axest)**2)==0) spgroup=9 97 98 case(8) 99 100 n_axest=(/0,0,0,0,0,0,1,1,0,0, 0,0,0,0,1,2,0,1,0,2, 0,0,0,0,0,0,0,0,0,0,0/) 101 if(sum((n_axes-n_axest)**2)==0) spgroup=36 102 103 n_axest=(/0,0,0,0,2,0,1,1,1,0, 0,0,0,0,1,1,0,0,0,1, 0,0,0,0,0,0,0,0,0,0,0/) 104 if(sum((n_axes-n_axest)**2)==0) spgroup=12 105 n_axest=(/0,0,0,0,2,0,1,1,1,0, 0,0,0,0,0,1,0,1,0,1, 0,0,0,0,0,0,0,0,0,0,0/) 106 if(sum((n_axes-n_axest)**2)==0) spgroup=15 107 n_axest=(/0,0,0,0,0,0,1,1,1,0, 0,0,0,0,2,1,0,1,0,1, 0,0,0,0,0,0,0,0,0,0,0/) 108 if(sum((n_axes-n_axest)**2)==0) spgroup=38 109 n_axest=(/0,0,0,0,0,0,1,1,1,0, 0,0,0,0,1,3,0,0,0,1, 0,0,0,0,0,0,0,0,0,0,0/) 110 if(sum((n_axes-n_axest)**2)==0) spgroup=39 111 n_axest=(/0,0,0,0,0,0,1,1,1,0, 0,0,0,0,1,1,0,2,0,1, 0,0,0,0,0,0,0,0,0,0,0/) 112 if(sum((n_axes-n_axest)**2)==0) spgroup=40 113 n_axest=(/0,0,0,0,0,0,1,1,1,0, 0,0,0,0,0,3,0,1,0,1, 0,0,0,0,0,0,0,0,0,0,0/) 114 if(sum((n_axes-n_axest)**2)==0) spgroup=41 115 116 n_axest=(/0,0,0,0,0,0,1,1,2,0, 0,0,0,0,0,0,0,0,0,4, 0,0,0,0,0,0,0,0,0,0,0/) 117 if(sum((n_axes-n_axest)**2)==0) spgroup=20 118 n_axest=(/0,0,0,0,0,0,1,1,2,0, 0,0,0,0,2,2,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0/) 119 if(sum((n_axes-n_axest)**2)==0) spgroup=35 120 n_axest=(/0,0,0,0,0,0,1,1,2,0, 0,0,0,0,0,2,0,2,0,0, 0,0,0,0,0,0,0,0,0,0,0/) 121 if(sum((n_axes-n_axest)**2)==0) spgroup=37 122 123 n_axest=(/0,0,0,0,0,0,1,1,4,0, 0,0,0,0,0,0,0,0,0,2, 0,0,0,0,0,0,0,0,0,0,0/) 124 if(sum((n_axes-n_axest)**2)==0) spgroup=21 125 126 case(16) 127 128 n_axest=(/0,0,0,0,2,0,1,1,2,0, 0,0,0,0,2,2,0,2,0,4, 0,0,0,0,0,0,0,0,0,0,0/) 129 if(sum((n_axes-n_axest)**2)==0) spgroup=63 130 n_axest=(/0,0,0,0,2,0,1,1,2,0, 0,0,0,0,1,4,0,1,0,4, 0,0,0,0,0,0,0,0,0,0,0/) 131 if(sum((n_axes-n_axest)**2)==0) spgroup=64 132 133 n_axest=(/0,0,0,0,2,0,1,1,4,0, 0,0,0,0,3,2,0,1,0,2, 0,0,0,0,0,0,0,0,0,0,0/) 134 if(sum((n_axes-n_axest)**2)==0) spgroup=65 135 n_axest=(/0,0,0,0,2,0,1,1,4,0, 0,0,0,0,2,4,0,0,0,2, 0,0,0,0,0,0,0,0,0,0,0/) 136 if(sum((n_axes-n_axest)**2)==0) spgroup=67 137 n_axest=(/0,0,0,0,2,0,1,1,4,0, 0,0,0,0,0,4,0,2,0,2, 0,0,0,0,0,0,0,0,0,0,0/) 138 if(sum((n_axes-n_axest)**2)==0) spgroup=68 139 n_axest=(/0,0,0,0,2,0,1,1,4,0, 0,0,0,0,1,2,0,3,0,2, 0,0,0,0,0,0,0,0,0,0,0/) 140 if(sum((n_axes-n_axest)**2)==0) spgroup=66 141 142 end select 143 144 else if(brvltt==7)then ! Rhombohedral lattice 145 146 select case(nsym) 147 148 case(3) 149 150 n_axest=(/0,0,0,0,0,0,0,1,0,2, 0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0/) 151 if(sum((n_axes-n_axest)**2)==0) spgroup=146 152 153 case(6) 154 155 n_axest=(/0,0,2,0,1,0,0,1,0,2, 0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0/) 156 if(sum((n_axes-n_axest)**2)==0) spgroup=148 157 n_axest=(/0,0,0,0,0,0,0,1,3,2, 0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0/) 158 if(sum((n_axes-n_axest)**2)==0) spgroup=155 159 n_axest=(/0,0,0,0,0,0,0,1,0,2, 0,0,0,0,3,0,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0/) 160 if(sum((n_axes-n_axest)**2)==0) spgroup=160 161 n_axest=(/0,0,0,0,0,0,0,1,0,2, 0,0,0,0,0,3,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0/) 162 if(sum((n_axes-n_axest)**2)==0) spgroup=161 163 164 case(12) 165 166 n_axest=(/0,0,2,0,1,0,0,1,3,2, 0,0,0,0,3,0,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0/) 167 if(sum((n_axes-n_axest)**2)==0) spgroup=166 168 n_axest=(/0,0,2,0,1,0,0,1,3,2, 0,0,0,0,0,3,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0/) 169 if(sum((n_axes-n_axest)**2)==0) spgroup=167 170 171 end select 172 173 end if 174 175 end subroutine symlist_others