TABLE OF CONTENTS


ABINIT/m_spgdata [ Modules ]

[ Top ] [ Modules ]

NAME

  m_spgdata

FUNCTION

COPYRIGHT

  Copyright (C) 2000-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_spgdata
28 
29  use defs_basis
30  use m_abicore
31  use m_errors
32 
33  use m_symtk,     only : symdet
34  use m_geometry,  only : xred2xcart
35 
36  implicit none
37 
38  private

m_spgdata/getptgroupma [ Functions ]

[ Top ] [ m_spgdata ] [ Functions ]

NAME

 getptgroupma

FUNCTION

 Return magnetic point group number from the full point group number
 and the point group number of the non-magnetic symmetry operations.
 The (normal) point group numbers are taken from
 The International Tables for Crystallography
 Volume A, 1983 Ed. Theo Hahn, D. Reidel Publishing Company
 The magnetic point group number are taken from
 The mathematical theory of symmetry in solids, Representation theory for point
 groups and space groups, 1972, C.J. Bradley and A.P.
 Cracknell, Clarendon Press, Oxford.
 In particular, see table 7.1 of the latter reference

INPUTS

 ptgroup = character(len=5) point group of all the symmetry operation
 ptgroupha = character(len=5) point group of the non-magnetic symmetry operation (halved point group)

OUTPUT

 ptgroupma = magnetic point group number

PARENTS

      symanal

CHILDREN

SOURCE

2329 subroutine getptgroupma(ptgroup,ptgroupha,ptgroupma)
2330 
2331 
2332 !This section has been created automatically by the script Abilint (TD).
2333 !Do not modify the following lines by hand.
2334 #undef ABI_FUNC
2335 #define ABI_FUNC 'getptgroupma'
2336 !End of the abilint section
2337 
2338  implicit none
2339 
2340 !Arguments ------------------------------------
2341 !scalars
2342  integer,intent(out) :: ptgroupma
2343  character(len=5),intent(in) :: ptgroup,ptgroupha
2344 
2345 ! *************************************************************************
2346 
2347 !DEBUG
2348 !write(std_out,*)' getptgroupma : enter '
2349 !write(std_out,*)' ptgroup="',ptgroup,'"'
2350 !write(std_out,*)' ptgroupha="',ptgroupha,'"'
2351 !ENDDEBUG
2352 
2353  ptgroupma=0
2354  select case (ptgroup)
2355  case("   -1")
2356    ptgroupma=1
2357  case("    2")
2358    ptgroupma=2
2359  case("   -2")
2360    ptgroupma=3
2361  case("  2/m")
2362    if(ptgroupha=="    2")ptgroupma=4
2363    if(ptgroupha=="   -2")ptgroupma=5
2364    if(ptgroupha=="   -1")ptgroupma=6
2365  case("  222")
2366    ptgroupma=7
2367  case("  mm2")
2368    if(ptgroupha=="    2")ptgroupma=8
2369    if(ptgroupha=="   -2")ptgroupma=9
2370  case("  mmm")
2371    if(ptgroupha=="  222")ptgroupma=10
2372    if(ptgroupha=="  mm2")ptgroupma=11
2373    if(ptgroupha=="  2/m")ptgroupma=12
2374  case("    4")
2375    ptgroupma=13
2376  case("   -4")
2377    ptgroupma=14
2378  case("  422")
2379    if(ptgroupha=="    4")ptgroupma=15
2380    if(ptgroupha=="  222")ptgroupma=16
2381  case("  4/m")
2382    if(ptgroupha=="    4")ptgroupma=17
2383    if(ptgroupha=="   -4")ptgroupma=18
2384    if(ptgroupha=="  2/m")ptgroupma=19
2385  case("  4mm")
2386    if(ptgroupha=="    4")ptgroupma=20
2387    if(ptgroupha=="  mm2")ptgroupma=21
2388  case(" -42m")
2389    if(ptgroupha=="   -4")ptgroupma=22
2390    if(ptgroupha=="  222")ptgroupma=23
2391    if(ptgroupha=="  mm2")ptgroupma=24
2392  case("4/mmm")
2393    if(ptgroupha=="  422")ptgroupma=25
2394    if(ptgroupha=="  4mm")ptgroupma=26
2395    if(ptgroupha=="  mmm")ptgroupma=27
2396    if(ptgroupha==" -42m")ptgroupma=28
2397    if(ptgroupha=="  4/m")ptgroupma=29
2398  case("   32")
2399    ptgroupma=30
2400  case("   3m")
2401    ptgroupma=31
2402  case("   -6")
2403    ptgroupma=32
2404  case(" -62m")
2405    if(ptgroupha=="   -6")ptgroupma=33
2406    if(ptgroupha=="   3m")ptgroupma=34
2407    if(ptgroupha=="   32")ptgroupma=35
2408  case("    6")
2409    ptgroupma=36
2410  case("   -3")
2411    ptgroupma=37
2412  case("  -3m")
2413    if(ptgroupha=="   -3")ptgroupma=38
2414    if(ptgroupha=="   3m")ptgroupma=39
2415    if(ptgroupha=="   32")ptgroupma=40
2416  case("  622")
2417    if(ptgroupha=="    6")ptgroupma=41
2418    if(ptgroupha=="   32")ptgroupma=42
2419  case("  6/m")
2420    if(ptgroupha=="    6")ptgroupma=43
2421    if(ptgroupha=="   -3")ptgroupma=44
2422    if(ptgroupha=="   -6")ptgroupma=45
2423  case("  6mm")
2424    if(ptgroupha=="    6")ptgroupma=46
2425    if(ptgroupha=="   3m")ptgroupma=47
2426  case("6/mmm")
2427    if(ptgroupha==" -62m")ptgroupma=48
2428    if(ptgroupha=="  -3m")ptgroupma=49
2429    if(ptgroupha=="  622")ptgroupma=50
2430    if(ptgroupha=="  6mm")ptgroupma=51
2431    if(ptgroupha=="  6/m")ptgroupma=52
2432  case("  m-3")
2433    ptgroupma=53
2434  case(" -43m")
2435    ptgroupma=54
2436  case("  432")
2437    ptgroupma=55
2438  case(" m-3m")
2439    if(ptgroupha=="  432")ptgroupma=56
2440    if(ptgroupha==" -43m")ptgroupma=57
2441    if(ptgroupha=="  m-3")ptgroupma=58
2442  end select
2443 
2444 !DEBUG
2445 !write(std_out,*)' getptgroupma : exit '
2446 !write(std_out,*)' ptgroupma="',ptgroupma,'"'
2447 !ENDDEBUG
2448 
2449 end subroutine getptgroupma

m_spgdata/prtspgroup [ Functions ]

[ Top ] [ m_spgdata ] [ Functions ]

NAME

 prtspgroup

FUNCTION

 Print the space group (first, the dataset)

INPUTS

  bravais(11)=characteristics of Bravais lattice (see symlatt.f)
  genafm(3)=generator of magnetic translations, in case of
            Shubnikov type IV magnetic groups (if zero, the group is
            not a type IV magnetic group)
  iout=unit number of output file
  jdtset= actual number of the dataset to be read
  ptgroupma=magnetic point group, in case of
            Shubnikov type III magnetic groups (if zero, the group is
            not a type III magnetic group)
  spgroup=space group number

OUTPUT

PARENTS

      memory_eval

CHILDREN

      ptgmadata,spgdata,wrtout,xred2xcart

SOURCE

 78 subroutine prtspgroup(bravais,genafm,iout,jdtset,ptgroupma,spgroup)
 79 
 80 
 81 !This section has been created automatically by the script Abilint (TD).
 82 !Do not modify the following lines by hand.
 83 #undef ABI_FUNC
 84 #define ABI_FUNC 'prtspgroup'
 85 !End of the abilint section
 86 
 87  implicit none
 88 
 89 !Arguments ------------------------------------
 90 !scalars
 91  integer,intent(in) :: iout,jdtset,ptgroupma,spgroup
 92 !arrays
 93  integer,intent(in) :: bravais(11)
 94  real(dp),intent(inout) :: genafm(3)
 95 
 96 !Local variables -------------------------------
 97 !scalars
 98  integer :: center,iholohedry,ii,shubnikov,spgaxor,spgorig,sporder,sumgen
 99  character(len=1) :: brvsb
100  character(len=10) :: ptgrpmasb
101  character(len=15) :: intsb,ptintsb,ptschsb,schsb
102  character(len=35) :: intsbl
103  character(len=500) :: message
104  character(len=80) :: bravais_name
105 !arrays
106  integer :: genafmint(3)
107  real(dp) :: genafmconv(3),rprimdconv(3,3)
108 
109 !*************************************************************************
110 
111 !DEBUG
112 !write(std_out,*)' prtspgroup : enter '
113 !write(std_out,*)' ptgroupma=',ptgroupma
114 !write(std_out,*)' genafm(:)=',genafm(:)
115 !ENDDEBUG
116 
117  center=bravais(2)
118  iholohedry=bravais(1)
119 
120 !Determine the magnetic type
121  shubnikov=1
122  if(ptgroupma/=0)shubnikov=3
123  if(sum(abs(genafm(:)))>tol6)then
124    shubnikov=4
125 !  Produce genafm in conventional axes,
126    rprimdconv(:,1)=bravais(3:5)
127    rprimdconv(:,2)=bravais(6:8)
128    rprimdconv(:,3)=bravais(9:11)
129    if(center/=0)rprimdconv(:,:)=rprimdconv(:,:)*half
130    call xred2xcart(1,rprimdconv,genafmconv,genafm)
131 !  Gives the associated translation, with components in the
132 !  interval ]-0.5,0.5] .
133    genafmconv(:)=genafmconv(:)-nint(genafmconv(:)-tol6)
134    do ii=1,3
135      genafmint(ii)=-1
136      if(abs(genafmconv(ii)-zero)<tol6)genafmint(ii)=0
137      if(abs(genafmconv(ii)-half)<tol6)genafmint(ii)=1
138    end do
139    if(minval(genafmint(:))==-1)then
140      write(message, '(3a,3es12.2,a)' )&
141 &     'The magnetic translation generator,',ch10,&
142 &     'genafmconv(:)=',genafmconv(:),&
143 &     'could not be identified.'
144      MSG_BUG(message)
145    end if
146  end if
147 
148 !Determine whether the space group can be printed
149  if(iholohedry<=0)then
150    if(jdtset/=0)then
151      write(message,'(a,a,i5,a)')ch10,&
152 &     ' DATASET',jdtset,' : the unit cell is not primitive'
153    else
154      write(message,'(a,a)')ch10,&
155 &     ' Symmetries : the unit cell is not primitive'
156    end if
157    call wrtout(std_out,message,'COLL')
158    call wrtout(iout,message,'COLL')
159  else if(spgroup==0)then
160    if(jdtset/=0)then
161      write(message,'(a,a,i5,a)')ch10,&
162 &     ' DATASET',jdtset,' : the space group has not been recognized'
163    else
164      write(message,'(a,a)')ch10,&
165 &     ' Symmetries : the space group has not been recognized'
166    end if
167    call wrtout(std_out,message,'COLL')
168    call wrtout(iout,message,'COLL')
169  else
170 
171 !  ------------------------------------------------------------------
172 !  The space group can be printed
173 
174 !  Determine the Bravais lattice
175 
176    bravais_name=' (the Bravais lattice could not be identified)'
177 
178    if(iholohedry==7)then ! Cubic
179 
180      if(center==0) then
181        if(shubnikov/=4)bravais_name='cP (primitive cubic)'
182        if(shubnikov==4)bravais_name='cP_I (primitive cubic, inner magnetic, #33)'
183      else if(center==-1) then
184        if(shubnikov/=4)bravais_name='cI (body-center cubic)' ! Only non-magnetic is possible
185      else if(center==-3) then
186        if(shubnikov/=4)bravais_name='cF (face-center cubic)'
187        if(shubnikov==4)bravais_name='cF_s (face-center cubic, simple cubic magnetic, #35)'
188      end if
189 
190    else if(iholohedry==4)then ! Tetragonal
191 
192      if(center==0) then
193        if(shubnikov/=4)bravais_name='tP (primitive tetrag.)'
194        if(shubnikov==4)then
195          sumgen=sum(genafmint(:))
196          if(sumgen==1)bravais_name='tP_c (primitive tetrag., c-magnetic, #23)'
197          if(sumgen==2)bravais_name='tP_C (primitive tetrag., C-magnetic, #24)'
198          if(sumgen==3)bravais_name='tP_I (primitive tetrag., centered magnetic, #25)'
199        end if
200      else if(center==-1)then
201        if(shubnikov/=4)bravais_name='tI (body-center tetrag.)'
202        if(shubnikov==4)bravais_name='tI_c (body-center tetrag., simple tetragonal magnetic, #27)'
203      end if
204 
205    else if(iholohedry==3)then ! Orthorhombic
206 
207      if(center==0) then
208        if(shubnikov/=4)bravais_name='oP (primitive ortho.)'
209        if(shubnikov==4)then
210          sumgen=sum(genafmint(:))
211          if(sumgen==1)then
212            if(genafmint(1)==1)bravais_name='oP_a (primitive ortho., a-magnetic, #11)'
213            if(genafmint(2)==1)bravais_name='oP_b (primitive ortho., b-magnetic, #11)'
214            if(genafmint(3)==1)bravais_name='oP_c (primitive ortho., c-magnetic, #11)'
215          else if(sumgen==2)then
216            if(genafmint(1)==0)bravais_name='oP_A (primitive ortho., A-magnetic, #12)'
217            if(genafmint(2)==0)bravais_name='oP_B (primitive ortho., B-magnetic, #12)'
218            if(genafmint(3)==0)bravais_name='oP_C (primitive ortho., C-magnetic, #12)'
219          else if(sumgen==3)then
220            bravais_name='oP_I (primitive ortho., centered magnetic, #13)'
221          end if
222        end if
223      else if(center==-1)then
224        if(shubnikov/=4)bravais_name='oI (body-center ortho.)'
225        if(shubnikov==4)bravais_name='oI_c (body-center ortho., simple ortho. magn., #21)'
226      else if(center==1 .or. center==2 .or. center==3)then
227        if(shubnikov/=4) bravais_name='oC (1-face-center ortho.)'
228        if(shubnikov==4)then
229          sumgen=sum(genafmint(:))
230          if(sumgen==1)then
231 !          One should work more to distinguish these magnetic groups
232            bravais_name='oC_(a,b,c) (1-face-cent. ortho., 1-magn., #15 or 16)'
233          else if(sumgen==2)then
234            bravais_name='oC_A (1-face-centered ortho., 1-face-magnetic, #17)'
235          else if(sumgen==3)then
236            bravais_name='oC_c (C-face-centered ortho., c-magnetic, #15)'
237          end if
238        end if
239      else if(center==-3)then
240        if(shubnikov/=4)bravais_name='oF (face-center ortho.)'
241        if(shubnikov==4)bravais_name='oF_s (face-center ortho., simple ortho. magnetic, #19)'
242      end if
243 
244    else if(iholohedry==6)then ! Hexagonal
245 
246      if(shubnikov/=4)bravais_name='hP (primitive hexag.)'
247      if(shubnikov==4)bravais_name='hP_c (primitive hexag., c-magnetic, #29)'
248 
249    else if(iholohedry==5)then ! Rhombohedral
250 
251      if(shubnikov/=4)bravais_name='hR (rhombohedral)'
252      if(shubnikov==4)bravais_name='hR_I (rhombohedral, centered magnetic, #31)'
253 
254    else if(iholohedry==2)then ! Monoclinic
255 
256      if(center==0)then
257        if(shubnikov/=4)bravais_name='mP (primitive monocl.)'
258        if(shubnikov==4)then
259          sumgen=sum(genafmint(:))
260          if(sumgen==1)then
261            if(genafmint(1)==1)bravais_name='mP_a (primitive monocl., a-magnetic, #5)'
262            if(genafmint(2)==1)bravais_name='mP_b (primitive monocl., b-magnetic, #4)'
263            if(genafmint(3)==1)bravais_name='mP_c (primitive monocl., c-magnetic, #5)'
264          else if(sumgen==2)then
265            if(genafmint(1)==0)bravais_name='mP_A (primitive monocl., A-magnetic, #6)'
266            if(genafmint(2)==0)bravais_name='mP_B (primitive monocl., B-magnetic, #6)'
267            if(genafmint(3)==0)bravais_name='mP_C (primitive monocl., C-magnetic, #6)'
268          end if
269        end if
270      else if(center==3)then
271        if(shubnikov/=4)bravais_name='mC (1-face-center monocl.)'
272        if(shubnikov==4)then
273          if(genafmint(3)==1)bravais_name='mC_c (C-face-center monocl., c-magnetic, #8)'
274          if(genafmint(3)/=1)bravais_name='mC_a (C-face-center monocl., a-magnetic, #9)'
275        end if
276      else if(center==-3)then
277        if(shubnikov/=4)bravais_name='(reduction of face-center)'
278      end if
279 
280    else if(iholohedry==1)then ! Triclinic
281 
282      if(shubnikov/=4)bravais_name='aP (primitive triclinic)'
283      if(shubnikov==4)bravais_name='aP_s (primitive triclinic, simple magnetic, #2)'
284 
285    end if
286 
287 !  Determine the symbol of the Fedorov space group
288    spgaxor=1 ; spgorig=1
289    call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig)
290 
291 !  Prepare print of the dataset, symmetry point group, Bravais lattice
292    if(shubnikov==1)then
293 
294      if(jdtset/=0)then
295        write(message,'(a,a,i5,a,a,a,a,i3,a,a,a)' )ch10,&
296 &       ' DATASET',jdtset,' : space group ',trim(brvsb),trim(intsb),' (#',spgroup,')',&
297 &       '; Bravais ',trim(bravais_name)
298      else
299        write(message,'(a,a,a,a,a,i3,a,a,a)' )ch10,&
300 &       ' Symmetries : space group ',trim(brvsb),trim(intsb),' (#',spgroup,')',&
301 &       '; Bravais ',trim(bravais_name)
302      end if
303      call wrtout(std_out,message,'COLL')
304      call wrtout(iout,message,'COLL')
305 
306    else if(shubnikov==3)then
307 
308      if(jdtset/=0)then
309        write(message,'(a,a,i5,a)' )ch10,&
310 &       ' DATASET',jdtset,' : magnetic group, Shubnikov type III '
311      else
312        write(message,'(2a)' )ch10,&
313 &       ' Magnetic group, Shubnikov type III '
314      end if
315      call wrtout(std_out,message,'COLL')
316      call wrtout(iout,message,'COLL')
317 
318      write(message,'(a,a,a,a,i3,a,a,a)' )&
319 &     ' Fedorov space group ',trim(brvsb),trim(intsb),' (#',spgroup,')',&
320 &     '; Bravais ',trim(bravais_name)
321      call wrtout(std_out,message,'COLL')
322      call wrtout(iout,message,'COLL')
323 
324      call ptgmadata(ptgroupma,ptgrpmasb)
325 
326      write(message,'(3a,i3,a)' )&
327 &     ' Magnetic point group ',trim(ptgrpmasb),' (#',ptgroupma,')'
328      call wrtout(std_out,message,'COLL')
329      call wrtout(iout,message,'COLL')
330 
331    else if(shubnikov==4)then
332 
333      if(jdtset/=0)then
334        write(message,'(a,a,i5,a)' )ch10,&
335 &       ' DATASET',jdtset,' : magnetic group, Shubnikov type IV '
336      else
337        write(message,'(2a)' )ch10,&
338 &       ' Magnetic group, Shubnikov type IV '
339      end if
340      call wrtout(std_out,message,'COLL')
341      call wrtout(iout,message,'COLL')
342 
343      write(message,'(a,a,a,a,i3,a)' )&
344 &     ' Fedorov space group ',trim(brvsb),trim(intsb),' (#',spgroup,')'
345      call wrtout(std_out,message,'COLL')
346      call wrtout(iout,message,'COLL')
347 
348      write(message,'(2a)' )&
349 &     ' Magnetic Bravais lattice ',trim(bravais_name)
350      call wrtout(std_out,message,'COLL')
351      call wrtout(iout,message,'COLL')
352 
353    end if
354  end if
355 
356 end subroutine prtspgroup

m_spgdata/ptgmadata [ Functions ]

[ Top ] [ m_spgdata ] [ Functions ]

NAME

 ptgmadata

FUNCTION

 Return magnetic point group symbol from the magnetic point group number
 The symbols and numbers are taken from  The Internationl Tables for Crystallography
 Volume A, 1983 Ed. Theo Hahn, D. Reidel Publishing Company and
 The mathematical theory of symmetry in solids, Representation theory for point
 groups and space groups, 1972, C.J. Bradley and A.P.
 Cracknell, Clarendon Press, Oxford.

INPUTS

 ptgroupma = space group number

OUTPUT

 ptgrpmasb= symbol

PARENTS

      prtspgroup

CHILDREN

SOURCE

2159 subroutine ptgmadata(ptgroupma,ptgrpmasb)
2160 
2161 
2162 !This section has been created automatically by the script Abilint (TD).
2163 !Do not modify the following lines by hand.
2164 #undef ABI_FUNC
2165 #define ABI_FUNC 'ptgmadata'
2166 !End of the abilint section
2167 
2168  implicit none
2169 
2170 !Arguments ------------------------------------
2171 !scalars
2172  integer,intent(in) :: ptgroupma
2173  character(len=10),intent(out) :: ptgrpmasb
2174 
2175 ! *************************************************************************
2176 
2177  select case (ptgroupma)
2178  case(1)
2179    ptgrpmasb="-1'"
2180  case(2)
2181    ptgrpmasb="2'"
2182  case(3)
2183    ptgrpmasb="m'"
2184  case(4)
2185    ptgrpmasb="2/m'"
2186  case(5)
2187    ptgrpmasb="2'/m"
2188  case(6)
2189    ptgrpmasb="2'/m'"
2190  case(7)
2191    ptgrpmasb="2'2'2"
2192  case(8)
2193    ptgrpmasb="m'm'2"
2194  case(9)
2195    ptgrpmasb="m'm2'"
2196  case(10)
2197    ptgrpmasb="m'm'm'"
2198  case(11)
2199    ptgrpmasb="mmm'"
2200  case(12)
2201    ptgrpmasb="m'm'm"
2202  case(13)
2203    ptgrpmasb="4'"
2204  case(14)
2205    ptgrpmasb="-4'"
2206  case(15)
2207    ptgrpmasb="42'2'"
2208  case(16)
2209    ptgrpmasb="4'22'"
2210  case(17)
2211    ptgrpmasb="4/m'"
2212  case(18)
2213    ptgrpmasb="4'/m'"
2214  case(19)
2215    ptgrpmasb="4'/m"
2216  case(20)
2217    ptgrpmasb="4m'm'"
2218  case(21)
2219    ptgrpmasb="4'mm'"
2220  case(22)
2221    ptgrpmasb="-42'm'"
2222  case(23)
2223    ptgrpmasb="-4'2m'"
2224  case(24)
2225    ptgrpmasb="-4'm2'"
2226  case(25)
2227    ptgrpmasb="4/m'm'm'"
2228  case(26)
2229    ptgrpmasb="4/m'mm"
2230  case(27)
2231    ptgrpmasb="4'/mmm'"
2232  case(28)
2233    ptgrpmasb="4'/m'm'm"
2234  case(29)
2235    ptgrpmasb="4/mm'm'"
2236  case(30)
2237    ptgrpmasb="32'"
2238  case(31)
2239    ptgrpmasb="3m'"
2240  case(32)
2241    ptgrpmasb="-6'"
2242  case(33)
2243    ptgrpmasb="-6m'2'"
2244  case(34)
2245    ptgrpmasb="-6'm2'"
2246  case(35)
2247    ptgrpmasb="-6'm'2"
2248  case(36)
2249    ptgrpmasb="6'"
2250  case(37)
2251    ptgrpmasb="-3'"
2252  case(38)
2253    ptgrpmasb="-3m'"
2254  case(39)
2255    ptgrpmasb="-3'm"
2256  case(40)
2257    ptgrpmasb="-3'm'"
2258  case(41)
2259    ptgrpmasb="62'2'"
2260  case(42)
2261    ptgrpmasb="6'2'2"
2262  case(43)
2263    ptgrpmasb="6/m'"
2264  case(44)
2265    ptgrpmasb="6'/m'"
2266  case(45)
2267    ptgrpmasb="6'/m"
2268  case(46)
2269    ptgrpmasb="6m'm'"
2270  case(47)
2271    ptgrpmasb="6'm'm"
2272  case(48)
2273    ptgrpmasb="6'/mmm'"
2274  case(49)
2275    ptgrpmasb="6'/m'm'm"
2276  case(50)
2277    ptgrpmasb="6/m'm'm'"
2278  case(51)
2279    ptgrpmasb="6/m'mm"
2280  case(52)
2281    ptgrpmasb="6/mm'm'"
2282  case(53)
2283    ptgrpmasb="m'3"
2284  case(54)
2285    ptgrpmasb="-4'3m'"
2286  case(55)
2287    ptgrpmasb="4'32'"
2288  case(56)
2289    ptgrpmasb="m'3m'"
2290  case(57)
2291    ptgrpmasb="m'3m"
2292  case(58)
2293    ptgrpmasb="m3m'"
2294  end select
2295 
2296 end subroutine ptgmadata

m_spgdata/spgdata [ Functions ]

[ Top ] [ m_spgdata ] [ Functions ]

NAME

 spgdata

FUNCTION

 Return point and space group data: Bravais lattice symbol,
 international symbol, Schonflies symbol, multiplicity
 The symbols are taken from  The International Tables for Crystallography
 Volume A, 1983 Ed. Theo Hahn, D. Reidel Publishing Company and
 The mathematical theory of symmetry in solids, Representation theory for point
 groups and space groups, 1972, C.J. Bradley and A.P.
 Cracknell, Clarendon Press, Oxford.

INPUTS

 spgroup = space group number
 spgorig = space group origin
 spgaxor = space group axis orientation

OUTPUT

 brvsb=Bravais lattice symbol (P, I, F, A, B, C, R)
 intsb=international symbol (like m3m, 222, 2_12_12_1)
 intsbl=international symbol in long format like P2_b = P121)
 ptintsb=International point group symbol
 ptschsb=Schoenflies point group symbol
 sporder=multiplicity of the space group
 schsb=Schoenflies symbol

NOTES

 brvsb, intsb, and schsb have been extensively checked, while
 more checking should be done for the others
 XG20160612 : in particular, at present it might be that spgaxor and spgorig are indetermined
 (e.g. spgaxor=-1;spgorig=-1) at input.
 When this has a bearing on some of the output variables (even brvsb or intsb !),
 these are mentioned as being X, unknown, or to be determined.

PARENTS

      m_ab7_symmetry,prt_cif,prtspgroup,symsgcube,symsghexa,symsgmono
      symsgortho,symsgtetra,symspgr

CHILDREN

SOURCE

 403 subroutine spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig)
 404 
 405 
 406 !This section has been created automatically by the script Abilint (TD).
 407 !Do not modify the following lines by hand.
 408 #undef ABI_FUNC
 409 #define ABI_FUNC 'spgdata'
 410 !End of the abilint section
 411 
 412  implicit none
 413 
 414 !Arguments ------------------------------------
 415 !scalars
 416  integer,intent(in) :: spgaxor,spgorig,spgroup
 417  integer,intent(out) :: sporder
 418  character(len=1),intent(out) :: brvsb
 419  character(len=15),intent(out) :: intsb,ptintsb,ptschsb,schsb
 420  character(len=35),intent(out) :: intsbl
 421 
 422 ! *************************************************************************
 423 
 424  intsbl="same"
 425 !defaults for case spgroup is not well defined (eg chkprim 0)
 426  brvsb="P"
 427  intsb="1"
 428  schsb="C1^1"
 429  sporder=1
 430 
 431  select case (spgroup)
 432  case(1)
 433    brvsb="P"; intsb="1"; schsb="C1^1"; sporder=1
 434  case(2)
 435    brvsb="P"; intsb="-1"; schsb="Ci^1"; sporder=2
 436  case(3)
 437    brvsb="P"; intsb="2"; schsb="C2^1"; sporder=2
 438    select case (spgaxor)
 439    case(1)
 440      intsbl="P 2 _b = P 1 2 1"
 441    case(2)
 442      intsbl="P 2_a = P 2 1 1"
 443    case(3)
 444      intsbl="P 2 _c = P 1 1 2"
 445    case default
 446      intsbl="intsbl to be determined"
 447    end select
 448  case(4)
 449    brvsb="P"; intsb="2_1"; schsb="C2^2"; sporder=2
 450    select case (spgaxor)
 451    case(1)
 452      intsbl="P 2 1 _b = P 1 2_1 1"
 453    case(2)
 454      intsbl="P 2 1 _a = P 2_1 1 1"
 455    case(3)
 456      intsbl="P 2 1 _c = P 1 1 2_1"
 457    case default
 458      intsbl="intsbl to be determined"
 459    end select
 460  case(5)
 461    brvsb="C"; intsb="2"; schsb="C2^3"; sporder=2
 462    select case (spgaxor)
 463    case(1)
 464      intsbl="C 2 _b1 =  C 1 2 1"
 465    case(2)
 466      intsbl="C 2 _a1 =  B 2 1 1"
 467    case(3)
 468      intsbl="C 2 _a2 =  C 2 1 1"
 469    case(4)
 470      intsbl="C 2 _a3 =  I 2 1 1"
 471    case(5)
 472      intsbl="C 2 _b2 =  A 1 2 1"
 473    case(6)
 474      intsbl="C 2 _b3 =  I 1 2 1"
 475    case(7)
 476      intsbl="C 2 _c1 =  A 1 1 2"
 477    case(8)
 478      intsbl="C 2 _c2 =  B 1 1 2 = B 2"
 479    case(9)
 480      intsbl="C 2 _c3 =  I 1 1 2"
 481    case default
 482      intsbl="intsbl to be determined"
 483    end select
 484  case(6)
 485    brvsb="P"; intsb="m"; schsb="Cs^1"; sporder=2
 486    select case (spgaxor)
 487    case(1)
 488      intsbl="P m _b = P 1 m 1"
 489    case(2)
 490      intsbl="P m _a = P m 1 1"
 491    case(3)
 492      intsbl="P m _c = P 1 1 m"
 493    case default
 494      intsbl="intsbl to be determined"
 495    end select
 496  case(7)
 497    brvsb="P"; intsb="c"; schsb="Cs^2"; sporder=2
 498    select case (spgaxor)
 499    case(1)
 500      intsbl="P c _b1 = P 1 c 1"
 501    case(2)
 502      intsbl="P c _a1 = P b 1 1"
 503    case(3)
 504      intsbl="P c _a2 = P n 1 1"
 505    case(4)
 506      intsbl="P c _a3 = P c 1 1"
 507    case(5)
 508      intsbl="P c _b2 = P 1 n 1"
 509    case(6)
 510      intsbl="P c _b3 = P 1 a 1"
 511    case(7)
 512      intsbl="P c _c1 = P 1 1 a"
 513    case(8)
 514      intsbl="P c _c2 = P 1 1 n"
 515    case(9)
 516      intsbl="P c _c3 = P 1 1 b = P b"
 517    case default
 518      intsbl="intsbl to be determined"
 519    end select
 520  case(8)
 521    brvsb="C"; intsb="m"; schsb="Cs^3"; sporder=4
 522    select case (spgaxor)
 523    case(1)
 524      intsbl="C m _b1 = C 1 m 1"
 525    case(2)
 526      intsbl="C m _a1 = B m 1 1"
 527    case(3)
 528      intsbl="C m _a2 = C m 1 1"
 529    case(4)
 530      intsbl="C m _a3 = I m 1 1"
 531    case(5)
 532      intsbl="C m _b2 = A 1 m 1"
 533    case(6)
 534      intsbl="C m _b3 = I 1 m 1"
 535    case(7)
 536      intsbl="C m _c1 = A 1 1 m"
 537    case(8)
 538      intsbl="C m _c2 = B 1 1 m = B m"
 539    case(9)
 540      intsbl="C m _c3 = I 1 1 m"
 541    case default
 542      intsbl="intsbl to be determined"
 543    end select
 544  case(9)
 545    brvsb="C"; intsb="c"; schsb="Cs^4"; sporder=4
 546    select case (spgaxor)
 547    case(1)
 548      intsbl="C c _b1 = C 1 c 1"
 549    case(2)
 550      intsbl="C c _a1 = B b 1 1"
 551    case(3)
 552      intsbl="C c _a2 = C n 1 1"
 553    case(4)
 554      intsbl="C c _a3 = I c 1 1"
 555    case(5)
 556      intsbl="C c _b2 = A 1 n 1"
 557    case(6)
 558      intsbl="C c _b3 = I 1 a 1"
 559    case(7)
 560      intsbl="C c _c1 = A 1 1 a"
 561    case(8)
 562      intsbl="C c _c2 = B 1 1 n"
 563    case(9)
 564      intsbl="C c _c3 = I 1 1 b"
 565    case default
 566      intsbl="intsbl to be determined"
 567    end select
 568  case(10)
 569    brvsb="P"; intsb="2/m"; schsb="C2h^1"; sporder=4
 570    select case (spgaxor)
 571    case(1)
 572      intsbl="P 2/m _b = P 1 2/m 1"
 573    case(2)
 574      intsbl="P 2/m _a = P 2/m 1 1"
 575    case(3)
 576      intsbl="P 2/m _c = P 1 1 2/m"
 577    case default
 578      intsbl="intsbl to be determined"
 579    end select
 580  case(11)
 581    brvsb="P"
 582    intsb="2_1/m"
 583    schsb="C2h^2"
 584    sporder=4
 585    select case (spgaxor)
 586    case(1)
 587      intsbl="P 2_1/m _b = P 1 2_1/m 1"
 588    case(2)
 589      intsbl="P 2_1/m _a = P 2_1/m 1 1"
 590    case(3)
 591      intsbl="P 2_1/m _c = P 1 1 2_1/m"
 592    case default
 593      intsbl="intsbl to be determined"
 594    end select
 595  case(12)
 596    brvsb="C"; intsb="2/m"; schsb="C2h^3"; sporder=8
 597    select case (spgaxor)
 598    case(1)
 599      intsbl="C 2/m _b1 = C 1 2/m 1"
 600    case(2)
 601      intsbl="C 2/m _a1 = B 2/m 1 1"
 602    case(3)
 603      intsbl="C 2/m _a2 = C 2/m 1 1"
 604    case(4)
 605      intsbl="C 2/m _a3 = I 2/m 1 1"
 606    case(5)
 607      intsbl="C 2/m _b2 = A 1 2/m 1"
 608    case(6)
 609      intsbl="C 2/m _b3 = I 1 2/m 1"
 610    case(7)
 611      intsbl="C 2/m _c1 = A 1 1 2/m"
 612    case(8)
 613      intsbl="C 2/m _c2 = B 1 1 2/m = B 2/m"
 614    case(9)
 615      intsbl="C 2/m _c3 = I 1 1 2/m"
 616    case default
 617      intsbl="intsbl to be determined"
 618    end select
 619  case(13)
 620    brvsb="P"; intsb="2/c"; schsb="C2h^4"; sporder=4
 621    select case (spgaxor)
 622    case(1)
 623      intsbl="P 2/c _b1 = P 1 2/c 1"
 624    case(2)
 625      intsbl="P 2/c _a1 = P 2/b 1 1"
 626    case(3)
 627      intsbl="P 2/c _a2 = P 2/n 1 1"
 628    case(4)
 629      intsbl="P 2/c _a3 = P 2/c 1 1"
 630    case(5)
 631      intsbl="P 2/c _b2 = P 1 2/n 1"
 632    case(6)
 633      intsbl="P 2/c _b3 = P 1 2/a 1"
 634    case(7)
 635      intsbl="P 2/c _c1 = P 1 1 2/a"
 636    case(8)
 637      intsbl="P 2/c _c2 = P 1 1 2/n"
 638    case(9)
 639      intsbl="P 2/c _c3 = P 1 1 2/b = P 2/b"
 640    case default
 641      intsbl="intsbl to be determined"
 642    end select
 643  case(14)
 644    brvsb="P"; intsb="2_1/c"; schsb="C2h^5"; sporder=4
 645    select case (spgaxor)
 646    case(1)
 647      intsbl="P 2_1/c _b1 = P 1 2_1/c 1"
 648    case(2)
 649      intsbl="P 2_1/c _a1 = P 2_1/b 1 1"
 650    case(3)
 651      intsbl="P 2_1/c _a2 = P 2_1/n 1 1"
 652    case(4)
 653      intsbl="P 2_1/c _a3 = P 2_1/c 1 1"
 654    case(5)
 655      intsbl="P 2_1/c _b2 = P 1 2_1/n 1"
 656    case(6)
 657      intsbl="P 2_1/c _b3 = P 1 2_1/a 1"
 658    case(7)
 659      intsbl="P 2_1/c _c1 = P 1 1 2_1/a"
 660    case(8)
 661      intsbl="P 2_1/c _c2 = P 1 1 2_1/n"
 662    case(9)
 663      intsbl="P 2_1/c _c3 = P 1 1 2_1/b = P 2_1/b"
 664    case default
 665      intsbl="intsbl to be determined"
 666    end select
 667  case(15)
 668    brvsb="C"; intsb="2/c"; schsb="C2h^6"; sporder=8
 669    select case (spgaxor)
 670    case(1)
 671      intsbl="C 2/c _b1 = C 1 2/c 1"
 672    case(2)
 673      intsbl="C 2/c _a1 = B 2/b 1 1"
 674    case(3)
 675      intsbl="C 2/c _a2 = C 2/n 1 1"
 676    case(4)
 677      intsbl="C 2/c _a3 = I 2/c 1 1"
 678    case(5)
 679      intsbl="C 2/c _b2 = A 1 2/n 1"
 680    case(6)
 681      intsbl="C 2/c _b3 = I 1 2/a 1"
 682    case(7)
 683      intsbl="C 2/c _c1 = A 1 1 2/a"
 684    case(8)
 685      intsbl="C 2/c _c2 = B 1 1 2/n"
 686    case(9)
 687      intsbl="C 2/c _c3 = I 1 1 2/b"
 688    case default
 689      intsbl="intsbl to be determined"
 690    end select
 691  case(16)
 692    brvsb="P"; intsb="2 2 2"; schsb="D2^1"; sporder=4
 693  case(17)
 694    brvsb="P"; intsb="2 2 2_1"; schsb="D2^2"; sporder=4
 695    select case (spgaxor)
 696    case(1)
 697      intsbl="P 2 2 2_1"
 698    case(2)
 699      intsbl="P 2_1 2 2"
 700    case(3)
 701      intsbl="P 2 2_1 2"
 702    case default
 703      intsbl="intsbl to be determined"
 704    end select
 705  case(18)
 706    brvsb="P"; intsb="2_1 2_1 2"; schsb="D2^3"; sporder=4
 707    select case (spgaxor)
 708    case(1)
 709      intsbl="P 2_1 2_1 2"
 710    case(2)
 711      intsbl="P 2 2_1 2_1"
 712    case(3)
 713      intsbl="P 2_1 2 2_1"
 714    case default
 715      intsbl="intsbl to be determined"
 716    end select
 717  case(19)
 718    brvsb="P"; intsb="2_1 2_1 2_1"; schsb="D2^4"; sporder=4
 719  case(20)
 720    schsb="D2^5"; sporder=8
 721    select case (spgaxor)
 722    case(1)
 723      brvsb="C"; intsb="2 2 2_1"
 724    case(2)
 725      brvsb="A"; intsb="2_1 2 2"
 726    case(3)
 727      brvsb="B"; intsb="2 2_1 2"
 728    case default
 729      brvsb="X"
 730      intsbl="intsbl to be determined"
 731    end select
 732  case(21)
 733    schsb="D2^6"; sporder=8
 734    select case (spgaxor)
 735    case(1)
 736      brvsb="C"; intsb="2 2 2"
 737    case(2)
 738      brvsb="A"; intsb="2 2 2"
 739    case(3)
 740      brvsb="B"; intsb="2 2 2"
 741    case default
 742      brvsb="X"
 743      intsbl="intsbl to be determined"
 744    end select
 745  case(22)
 746    brvsb="F"; intsb="2 2 2"; schsb="D2^7"; sporder=16
 747  case(23)
 748    brvsb="I"; intsb="2 2 2"; schsb="D2^8"; sporder=8
 749  case(24)
 750    brvsb="I"; intsb="2_1 2_1 2_1"; schsb="D2^9"; sporder=8
 751  case(25)
 752    brvsb="P"; schsb="C2v^1"; sporder=4
 753    select case (spgaxor)
 754    case(1)
 755      intsb="m m 2"
 756    case(2)
 757      intsb="2 m m"
 758    case(3)
 759      intsb="m 2 m"
 760    case default
 761      intsb="intsb unknown"
 762    end select
 763  case(26)
 764    brvsb="P"; schsb="C2v^2"; sporder=4
 765    select case (spgaxor)
 766    case(1)
 767      intsb="m c 2_1"
 768    case(2)
 769      intsb="2_1 m a"
 770    case(3)
 771      intsb="b 2_1 m"
 772    case(4)
 773      intsb="m 2_1 b"
 774    case(5)
 775      intsb="c m 2_1"
 776    case(6)
 777      intsb="2_1 a m"
 778    case default
 779      intsb="intsb unknown"
 780    end select
 781  case(27)
 782    brvsb="P"; schsb="C2v^3"; sporder=4
 783    select case (spgaxor)
 784    case(1)
 785      intsb="c c 2"
 786    case(2)
 787      intsb="2 a a"
 788    case(3)
 789      intsb="b 2 b"
 790    case default
 791      intsb="intsb unknown"
 792    end select
 793  case(28)
 794    brvsb="P"; schsb="C2v^4"; sporder=4
 795    select case (spgaxor)
 796    case(1)
 797      intsb="m a 2"
 798    case(2)
 799      intsb="2 m b"
 800    case(3)
 801      intsb="c 2 m"
 802    case(4)
 803      intsb="m 2 a"
 804    case(5)
 805      intsb="b m 2"
 806    case(6)
 807      intsb="2 c m"
 808    case default
 809      intsb="intsb unknown"
 810    end select
 811  case(29)
 812    brvsb="P"; schsb="C2v^5"; sporder=4
 813    select case (spgaxor)
 814    case(1)
 815      intsb="c a 2_1"
 816    case(2)
 817      intsb="2_1 a b"
 818    case(3)
 819      intsb="c 2_1 b"
 820    case(4)
 821      intsb="b 2_1 a"
 822    case(5)
 823      intsb="b c 2_1"
 824    case(6)
 825      intsb="2_1 c a"
 826    case default
 827      intsb="intsb unknown"
 828    end select
 829  case(30)
 830    brvsb="P"; schsb="C2v^6"; sporder=4
 831    select case (spgaxor)
 832    case(1)
 833      intsb="n c 2"
 834    case(2)
 835      intsb="2 n a"
 836    case(3)
 837      intsb="b 2 n"
 838    case(4)
 839      intsb="n 2 b"
 840    case(5)
 841      intsb="c n 2"
 842    case(6)
 843      intsb="2 a n"
 844    case default
 845      intsb="intsb unknown"
 846    end select
 847  case(31)
 848    brvsb="P"; schsb="C2v^7"; sporder=4
 849    select case (spgaxor)
 850    case(1)
 851      intsb="m n 2_1"
 852    case(2)
 853      intsb="2_1 m n"
 854    case(3)
 855      intsb="n 2_1 m"
 856    case(4)
 857      intsb="m 2_1 n"
 858    case(5)
 859      intsb="n m 2_1"
 860    case(6)
 861      intsb="2_1 n m"
 862    case default
 863      intsb="intsb unknown"
 864    end select
 865  case(32)
 866    brvsb="P"; schsb="C2v^8"; sporder=4
 867    select case (spgaxor)
 868    case(1)
 869      intsb="b a 2"
 870    case(2)
 871      intsb="2 c b"
 872    case(3)
 873      intsb="c 2 a"
 874    case default
 875      intsb="intsb unknown"
 876    end select
 877  case(33)
 878    brvsb="P"; schsb="C2v^9"; sporder=4
 879    select case (spgaxor)
 880    case(1)
 881      intsb="n a 2_1"
 882    case(2)
 883      intsb="2_1 n b"
 884    case(3)
 885      intsb="c 2_1 n"
 886    case(4)
 887      intsb="n 2_1 a"
 888    case(5)
 889      intsb="b n 2_1"
 890    case(6)
 891      intsb="2_1 c n"
 892    case default
 893      intsb="intsb unknown"
 894    end select
 895  case(34)
 896    brvsb="P"; schsb="C2v^10"; sporder=4
 897    select case (spgaxor)
 898    case(1)
 899      intsb="n n 2"
 900    case(2)
 901      intsb="2 n n"
 902    case(3)
 903      intsb="n 2 n"
 904    case default
 905      intsb="intsb unknown"
 906    end select
 907  case(35)
 908    schsb="C2v^11"; sporder=8
 909    select case (spgaxor)
 910    case(1)
 911      brvsb="C"; intsb="m m 2"
 912    case(2)
 913      brvsb="A"; intsb="2 m m"
 914    case(3)
 915      brvsb="B"; intsb="m 2 m"
 916    case default
 917      brvsb="X"
 918      intsb="intsb unknown"
 919    end select
 920  case(36)
 921    schsb="C2v^12"; sporder=8
 922    select case (spgaxor)
 923    case(1)
 924      brvsb="C"; intsb="m c 2_1"
 925    case(2)
 926      brvsb="A"; intsb="2_1 m a"
 927    case(3)
 928      brvsb="B"; intsb="b 2_1 m"
 929    case(4)
 930      brvsb="B"; intsb="m 2_1 b"
 931    case(5)
 932      brvsb="C"; intsb="c m 2_1"
 933    case(6)
 934      brvsb="A"; intsb="2_1 a m"
 935    case default
 936      brvsb="X"
 937      intsb="intsb unknown"
 938    end select
 939  case(37)
 940    schsb="C2v^13"; sporder=8
 941    select case (spgaxor)
 942    case(1)
 943      brvsb="C"; intsb="c c 2"
 944    case(2)
 945      brvsb="A"; intsb="2 a a"
 946    case(3)
 947      brvsb="B"; intsb="b 2 b"
 948    case default
 949      brvsb="X"
 950      intsb="intsb unknown"
 951    end select
 952  case(38)
 953    schsb="C2v^14"; sporder=8
 954    select case (spgaxor)
 955    case(1)
 956      brvsb="A"; intsb="m m 2"
 957    case(2)
 958      brvsb="B"; intsb="2 m m"
 959    case(3)
 960      brvsb="C"; intsb="m 2 m"
 961    case(4)
 962      brvsb="A"; intsb="m 2 m"
 963    case(5)
 964      brvsb="B"; intsb="m m 2"
 965    case(6)
 966      brvsb="C"; intsb="2 m m"
 967    case default
 968      brvsb="X"
 969      intsb="intsb unknown"
 970    end select
 971  case(39)
 972    schsb="C2v^15"; sporder=8
 973    select case (spgaxor)
 974    case(1)
 975      brvsb="A"; intsb="b m 2"
 976    case(2)
 977      brvsb="B"; intsb="2 c m"
 978    case(3)
 979      brvsb="C"; intsb="m 2 a"
 980    case(4)
 981      brvsb="A"; intsb="c 2 m"
 982    case(5)
 983      brvsb="B"; intsb="m a 2"
 984    case(6)
 985      brvsb="C"; intsb="2 m b"
 986    case default
 987      brvsb="X"
 988      intsb="intsb unknown"
 989    end select
 990  case(40)
 991    schsb="C2v^16"; sporder=8
 992    select case (spgaxor)
 993    case(1)
 994      brvsb="A"; intsb="m a 2"
 995    case(2)
 996      brvsb="B"; intsb="2 m b"
 997    case(3)
 998      brvsb="C"; intsb="c 2 m"
 999    case(4)
1000      brvsb="A"; intsb="m 2 a"
1001    case(5)
1002      brvsb="B"; intsb="b m 2"
1003    case(6)
1004      brvsb="C"; intsb="2 c m"
1005    case default
1006      brvsb="X"
1007      intsb="intsb unknown"
1008    end select
1009  case(41)
1010    schsb="C2v^17"; sporder=8
1011    select case (spgaxor)
1012    case(1)
1013      brvsb="A"; intsb="b a 2"
1014    case(2)
1015      brvsb="B"; intsb="2 c b"
1016    case(3)
1017      brvsb="C"; intsb="c 2 a"
1018    case(4)
1019      brvsb="A"; intsb="c 2 a"
1020    case(5)
1021      brvsb="B"; intsb="b a 2"
1022    case(6)
1023      brvsb="C"; intsb="2 c b"
1024    case default
1025      brvsb="X"
1026      intsb="intsb unknown"
1027    end select
1028  case(42)
1029    brvsb="F"; schsb="C2v^18"; sporder=16
1030    select case (spgaxor)
1031    case(1)
1032      intsb="m m 2"
1033    case(2)
1034      intsb="2 m m"
1035    case(3)
1036      intsb="m 2 m"
1037    case default
1038      intsb="intsb unknown"
1039    end select
1040  case(43)
1041    brvsb="F"; schsb="C2v^19"; sporder=16
1042    select case (spgaxor)
1043    case(1)
1044      intsb="d d 2"
1045    case(2)
1046      intsb="2 d d"
1047    case(3)
1048      intsb="d 2 d"
1049    case default
1050      intsb="intsb unknown"
1051    end select
1052  case(44)
1053    brvsb="I"; schsb="C2v^20"; sporder=8
1054    select case (spgaxor)
1055    case(1)
1056      intsb="m m 2"
1057    case(2)
1058      intsb="2 m m"
1059    case(3)
1060      intsb="m 2 m"
1061    case default
1062      intsb="intsb unknown"
1063    end select
1064  case(45)
1065    brvsb="I"; schsb="C2v^21"; sporder=8
1066    select case (spgaxor)
1067    case(1)
1068      intsb="b a 2"
1069    case(2)
1070      intsb="2 c b"
1071    case(3)
1072      intsb="c 2 a"
1073    case default
1074      intsb="intsb unknown"
1075    end select
1076  case(46)
1077    brvsb="I"; schsb="C2v^22"; sporder=8
1078    select case (spgaxor)
1079    case(1)
1080      intsb="m a 2"
1081    case(2)
1082      intsb="2 m b"
1083    case(3)
1084      intsb="c 2 m"
1085    case(4)
1086      intsb="m 2 a"
1087    case(5)
1088      intsb="b m 2"
1089    case(6)
1090      intsb="2 c m"
1091    case default
1092      intsb="intsb unknown"
1093    end select
1094  case(47)
1095    brvsb="P"; intsb="m m m"; schsb="D2h^1"; sporder=8
1096  case(48)
1097    brvsb="P"; intsb="n n n"; schsb="D2h^2"; sporder=8
1098    select case (spgorig)
1099    case(1)
1100      intsbl="n n n _1"
1101    case(2)
1102      intsbl="n n n _2"
1103    case default
1104      intsbl="intsbl to be determined"
1105    end select
1106  case(49)
1107    brvsb="P"; schsb="D2h^3"; sporder=8
1108    select case (spgaxor)
1109    case(1)
1110      intsb="c c m"
1111    case(2)
1112      intsb="m a a"
1113    case(3)
1114      intsb="b m b"
1115    case default
1116      intsb="intsb unknown"
1117    end select
1118  case(50)
1119    brvsb="P"; schsb="D2h^4"; sporder=8
1120    select case(spgorig)
1121    case(1)
1122      select case(spgaxor)
1123      case(1)
1124        intsb="b a n"; intsbl="b a n _1"
1125      case(2)
1126        intsb="n c b"; intsbl="n c b _1"
1127      case(3)
1128        intsb="c n a"; intsbl="c n a _1"
1129      case default
1130        intsb="intsb unknown"
1131        intsbl="intsbl to be determined"
1132      end select
1133    case(2)
1134      select case(spgaxor)
1135      case(5)
1136        intsb="b a n"; intsbl="b a n _2"
1137      case(6)
1138        intsb="n c b"; intsbl="n c b _2"
1139      case(4)
1140        intsb="c n a"; intsbl="c n a _2"
1141      case default
1142        intsb="intsb unknown"
1143        intsbl="intsbl to be determined"
1144      end select
1145    case default
1146      intsb="intsb unknown"
1147      intsbl="intsbl to be determined"
1148    end select
1149  case(51)
1150    brvsb="P"; schsb="D2h^5"; sporder=8
1151    select case (spgaxor)
1152    case(1)
1153      intsb="m m a"
1154    case(2)
1155      intsb="b m m"
1156    case(3)
1157      intsb="m c m"
1158    case(4)
1159      intsb="m a m"
1160    case(5)
1161      intsb="m m b"
1162    case(6)
1163      intsb="c m m"
1164    case default
1165      intsb="intsb unknown"
1166    end select
1167  case(52)
1168    brvsb="P"; schsb="D2h^6"; sporder=8
1169    select case (spgaxor)
1170    case(1)
1171      intsb="n n a"
1172    case(2)
1173      intsb="b n n"
1174    case(3)
1175      intsb="n c n"
1176    case(4)
1177      intsb="n a n"
1178    case(5)
1179      intsb="n n b"
1180    case(6)
1181      intsb="c n n"
1182    case default
1183      intsb="intsb unknown"
1184    end select
1185  case(53)
1186    brvsb="P"
1187    schsb="D2h^7"
1188    sporder=8
1189    select case (spgaxor)
1190    case(1)
1191      intsb="m n a"
1192    case(2)
1193      intsb="b m n"
1194    case(3)
1195      intsb="n c m"
1196    case(4)
1197      intsb="m a n"
1198    case(5)
1199      intsb="n m b"
1200    case(6)
1201      intsb="c n m"
1202    case default
1203      intsb="intsb unknown"
1204    end select
1205  case(54)
1206    brvsb="P"; schsb="D2h^8"; sporder=8
1207    select case (spgaxor)
1208    case(1)
1209      intsb="c c a"
1210    case(2)
1211      intsb="b a a"
1212    case(3)
1213      intsb="b c b"
1214    case(4)
1215      intsb="b a b"
1216    case(5)
1217      intsb="c c b"
1218    case(6)
1219      intsb="c a a"
1220    case default
1221      intsb="intsb unknown"
1222    end select
1223  case(55)
1224    brvsb="P"; schsb="D2h^9"; sporder=8
1225    select case (spgaxor)
1226    case(1)
1227      intsb="b a m"
1228    case(2)
1229      intsb="m c b"
1230    case(3)
1231      intsb="c m a"
1232    case default
1233      intsb="intsb unknown"
1234    end select
1235  case(56)
1236    brvsb="P"; schsb="D2h^10"; sporder=8
1237    select case (spgaxor)
1238    case(1)
1239      intsb="c c n"
1240    case(2)
1241      intsb="n a a"
1242    case(3)
1243      intsb="b n b"
1244    case default
1245      intsb="intsb unknown"
1246    end select
1247  case(57)
1248    brvsb="P"; schsb="D2h^11"; sporder=8
1249    select case (spgaxor)
1250    case(1)
1251      intsb="b c m"
1252    case(2)
1253      intsb="m c a"
1254    case(3)
1255      intsb="b m a"
1256    case(4)
1257      intsb="c m b"
1258    case(5)
1259      intsb="c a m"
1260    case(6)
1261      intsb="m a b"
1262    case default
1263      intsb="intsb unknown"
1264    end select
1265  case(58)
1266    brvsb="P"; schsb="D2h^12"; sporder=8
1267    select case (spgaxor)
1268    case(1)
1269      intsb="n n m"
1270    case(2)
1271      intsb="m n n"
1272    case(3)
1273      intsb="n m n"
1274    case default
1275      intsb="intsb unknown"
1276    end select
1277  case(59)
1278    brvsb="P"; schsb="D2h^13"; sporder=8
1279    if (spgorig==1) then
1280      select case (spgaxor)
1281      case(1)
1282        intsb="m m n"; intsbl="m m n _1"
1283      case(2)
1284        intsb="m m n"; intsbl="n m m _1"
1285      case(3)
1286        intsb="m m n"; intsbl="m n m _1"
1287      case default
1288        intsb="intsb unknown"
1289        intsbl="intsbl to be determined"
1290      end select
1291    else if(spgorig==2) then
1292      select case (spgaxor)
1293      case(5)
1294        intsb="m m n"; intsbl="m m n _2"
1295      case(6)
1296        intsb="m m n"; intsbl="n m m _2"
1297      case(4)
1298        intsb="m m n"; intsbl="m n m _2"
1299      case default
1300        intsb="intsb unknown"
1301        intsbl="intsbl to be determined"
1302      end select
1303    else
1304      intsb="intsb unknown"
1305      intsbl="intsbl to be determined"
1306    end if
1307  case(60)
1308    brvsb="P"; schsb="D2h^14"; sporder=8
1309    select case (spgaxor)
1310    case(1)
1311      intsb="b c n"
1312    case(2)
1313      intsb="n c a"
1314    case(3)
1315      intsb="b n a"
1316    case(4)
1317      intsb="c n b"
1318    case(5)
1319      intsb="c a n"
1320    case(6)
1321      intsb="n a b"
1322    case default
1323      intsb="intsb unknown"
1324    end select
1325  case(61)
1326    brvsb="P"; schsb="D2h^15"; sporder=8
1327    if (spgaxor==1)then
1328      intsb="b c a"
1329    else if (spgaxor==2)then
1330      intsb="c a b"
1331    else
1332      intsb="intsb unknown"
1333    end if
1334  case(62)
1335    brvsb="P"; schsb="D2h^16"; sporder=8
1336    select case (spgaxor)
1337    case(1)
1338      intsb="n m a"
1339    case(2)
1340      intsb="b n m"
1341    case(3)
1342      intsb="m c n"
1343    case(4)
1344      intsb="n a m"
1345    case(5)
1346      intsb="m n b"
1347    case(6)
1348      intsb="c m n"
1349    case default
1350      intsb="intsb unknown"
1351    end select
1352  case(63)
1353    schsb="D2h^17"; sporder=16
1354    select case (spgaxor)
1355    case(1)
1356      brvsb="C"; intsb="m c m"
1357    case(2)
1358      brvsb="A"; intsb="m m a"
1359    case(3)
1360      brvsb="B"; intsb="b m m"
1361    case(4)
1362      brvsb="B"; intsb="m m b"
1363    case(5)
1364      brvsb="C"; intsb="c m m"
1365    case(6)
1366      brvsb="A"; intsb="m a m"
1367    case default
1368      brvsb="X"
1369      intsbl="intsbl unknown"
1370    end select
1371  case(64)
1372    schsb="D2h^18"; sporder=16
1373    select case (spgaxor)
1374    case(1)
1375      brvsb="C"; intsb="m c a"
1376    case(2)
1377      brvsb="A"; intsb="b m a"
1378    case(3)
1379      brvsb="B"; intsb="b c m"
1380    case(4)
1381      brvsb="B"; intsb="m a b"
1382    case(5)
1383      brvsb="C"; intsb="c m b"
1384    case(6)
1385      brvsb="A"; intsb="c a m"
1386    case default
1387      brvsb="X"
1388      intsb="intsb unknown"
1389    end select
1390  case(65)
1391    schsb="D2h^19"; sporder=16
1392    select case (spgaxor)
1393    case(1)
1394      brvsb="C"; intsb="m m m"
1395    case(2)
1396      brvsb="A"; intsb="m m m"
1397    case(3)
1398      brvsb="B"; intsb="m m m"
1399    case default
1400      brvsb="X"
1401      intsb="intsb unknown"
1402    end select
1403  case(66)
1404    schsb="D2h^20"; sporder=16
1405    select case (spgaxor)
1406    case(1)
1407      brvsb="C"; intsb="c c m"
1408    case(2)
1409      brvsb="A"; intsb="m a a"
1410    case(3)
1411      brvsb="B"; intsb="b m b"
1412    case default
1413      brvsb="X"
1414      intsb="intsb unknown"
1415    end select
1416  case(67)
1417    schsb="D2h^21"; sporder=16
1418    select case (spgaxor)
1419    case(1)
1420      brvsb="C"; intsb="m m a"
1421    case(2)
1422      brvsb="A"; intsb="b m m"
1423    case(3)
1424      brvsb="B"; intsb="m c m"
1425    case(4)
1426      brvsb="B"; intsb="m a m"
1427    case(5)
1428      brvsb="C"; intsb="m m b"
1429    case(6)
1430      brvsb="A"; intsb="c m m"
1431    case default
1432      brvsb="X"
1433      intsb="intsb unknown"
1434    end select
1435  case(68)
1436    schsb="D2h^22"; sporder=16
1437    if (spgorig==1) then
1438      select case (spgaxor)
1439      case(1)
1440        brvsb="C"; intsb="c c a"; intsbl="c c a _1"
1441      case(2)
1442        brvsb="A"; intsb="b a a"; intsbl="b a a _1"
1443      case(3)
1444        brvsb="B"; intsb="b c b"; intsbl="b c b _1"
1445      case(4)
1446        brvsb="B"; intsb="b a b"; intsbl="b a b _1"
1447      case(5)
1448        brvsb="C"; intsb="c c b"; intsbl="c c b _1"
1449      case(6)
1450        brvsb="A"; intsb="c a a"; intsbl="c a a _1"
1451      case default
1452        brvsb="X"
1453        intsb="intsb unknown"
1454        intsbl="intsbl to be determined"
1455      end select
1456    else if(spgorig==2)then
1457      select case (spgaxor)
1458      case(1)
1459        brvsb="C"; intsb="c c a"; intsbl="c c a _2"
1460      case(2)
1461        brvsb="A"; intsb="b a a"; intsbl="b a a _2"
1462      case(3)
1463        brvsb="B"; intsb="b c b"; intsbl="b c b _2"
1464      case(4)
1465        brvsb="B"; intsb="b a b"; intsbl="b a b _2"
1466      case(5)
1467        brvsb="C"; intsb="c c b"; intsbl="c c b _2"
1468      case(6)
1469        brvsb="A"; intsb="c a a"; intsbl="c a a _2"
1470      case default
1471        brvsb="X"
1472        intsb="intsb unknown"
1473        intsbl="intsbl to be determined"
1474      end select
1475    else
1476      brvsb="X"
1477      intsb="intsb unknown"
1478      intsbl="intsbl to be determined"
1479    end if
1480  case(69)
1481    brvsb="F"; intsb="m m m"; schsb="D2h^23"; sporder=32
1482  case(70)
1483    brvsb="F"; intsb="d d d"; schsb="D2h^24"; sporder=32
1484    if (spgorig==1)then
1485      intsbl="d d d _1"
1486    else if (spgorig==2)then
1487      intsbl="d d d _2"
1488    else
1489      intsbl="intsbl to be determined"
1490    end if
1491  case(71)
1492    brvsb="I"; intsb="m m m"; schsb="D2h^25"; sporder=16
1493  case(72)
1494    brvsb="I"; schsb="D2h^26"; sporder=16
1495    select case (spgaxor)
1496    case(1)
1497      intsb="b a m"
1498    case(2)
1499      intsb="m c b"
1500    case(3)
1501      intsb="c m a"
1502    case default
1503      intsb="intsb unknown"
1504    end select
1505  case(73)
1506    brvsb="I"; schsb="D2h^27"; sporder=16
1507    if (spgorig==1)then
1508      intsb="b c a"
1509    else if (spgorig==2)then
1510      intsb="c a b"
1511    else
1512      intsb="intsb unknown"
1513    end if
1514  case(74)
1515    brvsb="I"; schsb="D2h^28"; sporder=16
1516    select case (spgaxor)
1517    case(1)
1518      intsb="m m a"
1519    case(2)
1520      intsb="b m m"
1521    case(3)
1522      intsb="m c m"
1523    case(4)
1524      intsb="m a m"
1525    case(5)
1526      intsb="m m b"
1527    case(6)
1528      intsb="c m m"
1529    case default
1530      intsb="intsb unknown"
1531    end select
1532  case(75)
1533    brvsb="P"; intsb="4"; schsb="C4^1"; sporder=4
1534  case(76)
1535    brvsb="P"; intsb="4_1"; schsb="C4^2"; sporder=4
1536  case(77)
1537    brvsb="P"; intsb="4_2"; schsb="C4^3"; sporder=4
1538  case(78)
1539    brvsb="P"; intsb="4_3"; schsb="C4^4"; sporder=4
1540  case(79)
1541    brvsb="I"; intsb="4"; schsb="C4^5"; sporder=8
1542  case(80)
1543    brvsb="I"; intsb="4_1"; schsb="C4^6"; sporder=8
1544  case(81)
1545    brvsb="P"; intsb="-4"; schsb="S4^1"; sporder=4
1546  case(82)
1547    brvsb="I"; intsb="-4"; schsb="S4^2"; sporder=8
1548  case(83)
1549    brvsb="P"; intsb="4/m"; schsb="C4h^1"; sporder=8
1550  case(84)
1551    brvsb="P"; intsb="4_2/m"; schsb="C4h^2"; sporder=8
1552  case(85)
1553    brvsb="P"; intsb="4/n"; schsb="C4h^3"; sporder=8
1554    if (spgorig==1)then
1555      intsbl="4/n _1"
1556    else if (spgorig==2)then
1557      intsbl="4/n _2"
1558    else
1559      intsbl="intsbl to be determined"
1560    end if
1561  case(86)
1562    brvsb="P"; intsb="4_2/n"; schsb="C4h^4"; sporder=8
1563    if (spgorig==1)then
1564      intsbl="4_2/n _1"
1565    else if (spgorig==2)then
1566      intsbl="4_2/n _2"
1567    else
1568      intsbl="intsbl to be determined"
1569    end if
1570  case(87)
1571    brvsb="I"; intsb="4/m"; schsb="C4h^5"; sporder=16
1572  case(88)
1573    brvsb="I"; intsb="4_1/a"; schsb="C4h^6"; sporder=16
1574    if (spgorig==1)then
1575      intsbl="4_1/a _1"
1576    else if (spgorig==2)then
1577      intsbl="4_1/a _2"
1578    else
1579      intsbl="intsbl to be determined"
1580    end if
1581  case(89)
1582    brvsb="P"; intsb="4 2 2"; schsb="D4^1"; sporder=8
1583  case(90)
1584    brvsb="P"; intsb="4 2_1 2"; schsb="D4^2"; sporder=8
1585  case(91)
1586    brvsb="P"; intsb="4_1 2 2"; schsb="D4^3"; sporder=8
1587  case(92)
1588    brvsb="P"; intsb="4_1 2_1 2"; schsb="D4^4"; sporder=8
1589  case(93)
1590    brvsb="P"; intsb="4_2 2 2"; schsb="D4^5"; sporder=8
1591  case(94)
1592    brvsb="P"; intsb="4_2 2_1 2"; schsb="D4^6"; sporder=8
1593  case(95)
1594    brvsb="P"; intsb="4_3 2 2"; schsb="D4^7"; sporder=8
1595  case(96)
1596    brvsb="P"; intsb="4_3 2_1 2"; schsb="D4^8"; sporder=8
1597  case(97)
1598    brvsb="I"; intsb="4 2 2"; schsb="D4^9"; sporder=16
1599  case(98)
1600    brvsb="I"; intsb="4_1 2 2"; schsb="D4^10"; sporder=16
1601  case(99)
1602    brvsb="P"; intsb="4 m m"; schsb="C4v^1"; sporder=8
1603  case(100)
1604    brvsb="P"; intsb="4 b m"; schsb="C4v^2"; sporder=8
1605  case(101)
1606    brvsb="P"; intsb="4_2 c m"; schsb="C4v^3"; sporder=8
1607  case(102)
1608    brvsb="P"; intsb="4_2 n m"; schsb="C4v^4"; sporder=8
1609  case(103)
1610    brvsb="P"; intsb="4 c c"; schsb="C4v^5"; sporder=8
1611  case(104)
1612    brvsb="P"; intsb="4 n c"; schsb="C4v^6"; sporder=8
1613  case(105)
1614    brvsb="P"; intsb="4_2 m c"; schsb="C4v^7"; sporder=8
1615  case(106)
1616    brvsb="P"; intsb="4_2 b c"; schsb="C4v^8"; sporder=8
1617  case(107)
1618    brvsb="I"; intsb="4 m m"; schsb="C4v^9"; sporder=16
1619  case(108)
1620    brvsb="I"; intsb="4 c m"; schsb="C4v^10"; sporder=16
1621  case(109)
1622    brvsb="I"; intsb="4_1 m d"; schsb="C4v^11"; sporder=16
1623  case(110)
1624    brvsb="I"; intsb="4_1 c d"; schsb="C4v^12"; sporder=16
1625  case(111)
1626    brvsb="P"; intsb="-4 2 m"; schsb="D2d^1"; sporder=8
1627  case(112)
1628    brvsb="P"; intsb="-4 2 c"; schsb="D2d^2"; sporder=8
1629  case(113)
1630    brvsb="P"; intsb="-4 2_1 m"; schsb="D2d^3"; sporder=8
1631  case(114)
1632    brvsb="P"; intsb="-4 2_1 c"; schsb="D2d^4"; sporder=8
1633  case(115)
1634    brvsb="P"; intsb="-4 m 2"; schsb="D2d^5"; sporder=8
1635  case(116)
1636    brvsb="P"; intsb="-4 c 2"; schsb="D2d^6"; sporder=8
1637  case(117)
1638    brvsb="P"; intsb="-4 b 2"; schsb="D2d^7"; sporder=8
1639  case(118)
1640    brvsb="P"; intsb="-4 n 2"; schsb="D2d^8"; sporder=8
1641  case(119)
1642    brvsb="I"; intsb="-4 m 2"; schsb="D2d^9"; sporder=16
1643  case(120)
1644    brvsb="I"; intsb="-4 c 2"; schsb="D2d^10"; sporder=16
1645  case(121)
1646    brvsb="I"; intsb="-4 2 m"; schsb="D2d^11"; sporder=16
1647  case(122)
1648    brvsb="I"; intsb="-4 2 d"; schsb="D2d^12"; sporder=16
1649  case(123)
1650    brvsb="P"; intsb="4/m m m"; schsb="D4h^1"; sporder=16
1651  case(124)
1652    brvsb="P"; intsb="4/m c c"; schsb="D4h^2"; sporder=16
1653  case(125)
1654    brvsb="P"; intsb="4/n b m"; schsb="D4h^3"; sporder=16
1655    if (spgorig==1)then
1656      intsbl="4/n b m _1"
1657    else if (spgorig==2)then
1658      intsbl="4/n b m _2"
1659    else
1660      intsbl="intsbl to be determined"
1661    end if
1662  case(126)
1663    brvsb="P"; intsb="4/n n c"; schsb="D4h^4"; sporder=16
1664    if (spgorig==1)then
1665      intsbl="4/n n c _1"
1666    else if (spgorig==2)then
1667      intsbl="4/n n c _2"
1668    else
1669      intsbl="intsbl to be determined"
1670    end if
1671  case(127)
1672    brvsb="P"; intsb="4/m b m"; schsb="D4h^5"; sporder=16
1673  case(128)
1674    brvsb="P"; intsb="4/m n c"; schsb="D4h^6"; sporder=16
1675  case(129)
1676    brvsb="P"; intsb="4/n m m"; schsb="D4h^7"; sporder=16
1677    if (spgorig==1)then
1678      intsbl="4/n m m _1"
1679    else if (spgorig==2)then
1680      intsbl="4/n m m _2"
1681    else
1682      intsbl="intsbl to be determined"
1683    end if
1684  case(130)
1685    brvsb="P"; intsb="4/n c c"; schsb="D4h^8"; sporder=16
1686    if (spgorig==1)then
1687      intsbl="4/n c c _1"
1688    else if (spgorig==2) then
1689      intsbl="4/n c c _2"
1690    else
1691      intsbl="intsbl to be determined"
1692    end if
1693  case(131)
1694    brvsb="P"; intsb="4_2/m m c"; schsb="D4h^9"; sporder=16
1695  case(132)
1696    brvsb="P"; intsb="4_2/m c m"; schsb="D4h^10"; sporder=16
1697  case(133)
1698    brvsb="P"; intsb="4_2/n b c"; schsb="D4h^11"; sporder=16
1699    if (spgorig==1)then
1700      intsbl="4_2/n b c _1"
1701    else if (spgorig==2)then
1702      intsbl="4_2/n b c _2"
1703    else
1704      intsbl="intsbl to be determined"
1705    end if
1706  case(134)
1707    brvsb="P"; intsb="4_2/n n m"; schsb="D4h^12"; sporder=16
1708    if (spgorig==1)then
1709      intsbl="4_2/n n m _1"
1710    else if (spgorig==2)then
1711      intsbl="4_2/n n m _2"
1712    else
1713      intsbl="intsbl to be determined"
1714    end if
1715  case(135)
1716    brvsb="P"; intsb="4_2/m b c"; schsb="D4h^13"; sporder=16
1717  case(136)
1718    brvsb="P"; intsb="4_2/m n m"; schsb="D4h^14"; sporder=16
1719  case(137)
1720    brvsb="P"; intsb="4_2/n m c"; schsb="D4h^15"; sporder=16
1721    if (spgorig==1)then
1722      intsbl="4_2/n m c _1"
1723    else if (spgorig==2)then
1724      intsbl="4_2/n m c _2"
1725    else
1726      intsbl="intsbl to be determined"
1727    end if
1728  case(138)
1729    brvsb="P"; intsb="4_2/n c m"; schsb="D4h^16"; sporder=16
1730    if (spgorig==1)then
1731      intsbl="4_2/n c m _1"
1732    else if (spgorig==2)then
1733      intsbl="4_2/n c m _2"
1734    else
1735      intsbl="intsbl to be determined"
1736    end if
1737  case(139)
1738    brvsb="I"; intsb="4/m m m"; schsb="D4h^17"; sporder=32
1739  case(140)
1740    brvsb="I"; intsb="4/m c m"; schsb="D4h^18"; sporder=32
1741  case(141)
1742    brvsb="I"; intsb="4_1/a m d"; schsb="D4h^19"; sporder=32
1743    if (spgorig==1)then
1744      intsbl="4_1/a m d _1"
1745    else if (spgorig==2)then
1746      intsbl="4_1/a m d _2"
1747    else
1748      intsbl="intsbl to be determined"
1749    end if
1750  case(142)
1751    brvsb="I"; intsb="4_1/a c d"; schsb="D4h^20"; sporder=32
1752    if (spgorig==1)then
1753      intsbl="4_1/a c d _1"
1754    else if (spgorig==2)then
1755      intsbl="4_1/a c d _2"
1756    else
1757      intsbl="intsbl to be determined"
1758    end if
1759  case(143)
1760    brvsb="P"; intsb="3"; schsb="C3^1"; sporder=3
1761  case(144)
1762    brvsb="P"; intsb="3_1"; schsb="C3^2"; sporder=3
1763  case(145)
1764    brvsb="P"; intsb="3_2"; schsb="C3^3"; sporder=3
1765  case(146)
1766    brvsb="R"; intsb="3"; schsb="C3^4"
1767    if (spgorig==1)then
1768      intsbl="3 _H" ; sporder=9
1769    else if (spgorig==2)then
1770      intsbl="3 _R" ; sporder=3
1771    else
1772      intsbl="intsbl to be determined"
1773    end if
1774  case(147)
1775    brvsb="P"; intsb="-3"; schsb="C3i^1"; sporder=6
1776  case(148)
1777    brvsb="R"; intsb="-3"; schsb="C3i^2"
1778    if (spgorig==1) then
1779      intsbl="-3 _H" ; sporder=9
1780    else if (spgorig==2) then
1781      intsbl="-3 _R" ; sporder=3
1782    else
1783      intsbl="intsbl to be determined"
1784    end if
1785  case(149)
1786    brvsb="P"; intsb="3 1 2"; schsb="D3^1"; sporder=6
1787  case(150)
1788    brvsb="P"; intsb="3 2 1"; schsb="D3^2"; sporder=6
1789  case(151)
1790    brvsb="P"; intsb="3_1 1 2"; schsb="D3^3"; sporder=6
1791  case(152)
1792    brvsb="P"; intsb="3_1 2 1"; schsb="D3^4"; sporder=6
1793  case(153)
1794    brvsb="P"; intsb="3_2 1 2"; schsb="D3^5"; sporder=6
1795  case(154)
1796    brvsb="P"; intsb="3_2 2 1"; schsb="D3^6"; sporder=6
1797  case(155)
1798    brvsb="R"; intsb="3 2"; schsb="D3^7"
1799    if (spgorig==1) then
1800      intsbl="3 2 _H" ; sporder=18
1801    else if (spgorig==2) then
1802      intsbl="3 2 _R" ; sporder=6
1803    else
1804      intsbl="intsbl to be determined"
1805    end if
1806  case(156)
1807    brvsb="P"; intsb="3 m 1"; schsb="C3v^1"; sporder=6
1808  case(157)
1809    brvsb="P"; intsb="3 1 m"; schsb="C3v^2"; sporder=6
1810  case(158)
1811    brvsb="P"; intsb="3 c 1"; schsb="C3v^3"; sporder=6
1812  case(159)
1813    brvsb="P"; intsb="3 1 c"; schsb="C3v^4"; sporder=6
1814  case(160)
1815    brvsb="R"; intsb="3 m"; schsb="C3v^5"
1816    if (spgorig==1) then
1817      intsbl="3 m _H" ; sporder=18
1818    else if (spgorig==2) then
1819      intsbl="3 m _R" ; sporder=6
1820    else
1821      intsbl="intsbl to be determined"
1822    end if
1823  case(161)
1824    brvsb="R"; intsb="3 c"; schsb="C3v^6"
1825    if (spgorig==1) then
1826      intsbl="3 m _H" ; sporder=18
1827    else if (spgorig==2)then
1828      intsbl="3 m _R" ; sporder=6
1829    else
1830      intsbl="intsbl to be determined"
1831    end if
1832  case(162)
1833    brvsb="P"; intsb="-3 1 m"; schsb="D3d^1"; sporder=12
1834  case(163)
1835    brvsb="P"; intsb="-3 1 c"; schsb="D3d^2"; sporder=12
1836  case(164)
1837    brvsb="P"; intsb="-3 m 1"; schsb="D3d^3"; sporder=12
1838  case(165)
1839    brvsb="P"; intsb="-3 c 1"; schsb="D3d^4"; sporder=12
1840  case(166)
1841    brvsb="R"; intsb="-3 m"; schsb="D3d^5"
1842    if (spgorig==1) then
1843      intsbl="3 m _H"; sporder=18
1844    else if (spgorig==2) then
1845      intsbl="3 m _R"; sporder=6
1846    else
1847      intsbl="intsbl to be determined"
1848    end if
1849  case(167)
1850    brvsb="R"; intsb="-3 c"; schsb="D3d^6"
1851    if (spgorig==1) then
1852      intsbl="-3 c _H"; sporder=36
1853    else if (spgorig==2) then
1854      intsbl="-3 c _R"; sporder=12
1855    else
1856      intsbl="intsbl to be determined"
1857      sporder=-1
1858    end if
1859  case(168)
1860    brvsb="P"; intsb="6"; schsb="C6^1"; sporder=6
1861  case(169)
1862    brvsb="P"; intsb="6_1"; schsb="C6^2"; sporder=6
1863  case(170)
1864    brvsb="P"; intsb="6_5"; schsb="C6^3"; sporder=6
1865  case(171)
1866    brvsb="P"; intsb="6_2"; schsb="C6^4"; sporder=6
1867  case(172)
1868    brvsb="P"; intsb="6_4"; schsb="C6^5"; sporder=6
1869  case(173)
1870    brvsb="P"; intsb="6_3"; schsb="C6^6"; sporder=6
1871  case(174)
1872    brvsb="P"; intsb="-6"; schsb="C3h^1"; sporder=6
1873  case(175)
1874    brvsb="P"; intsb="6/m"; schsb="C6h^1"; sporder=12
1875  case(176)
1876    brvsb="P"; intsb="6_3/m"; schsb="C6h^2"; sporder=12
1877  case(177)
1878    brvsb="P"; intsb="6 2 2"; schsb="D6^1"; sporder=12
1879  case(178)
1880    brvsb="P"; intsb="6_1 2 2"; schsb="D6^2"; sporder=12
1881  case(179)
1882    brvsb="P"; intsb="6_5 2 2"; schsb="D6^3"; sporder=12
1883  case(180)
1884    brvsb="P"; intsb="6_2 2 2"; schsb="D6^4"; sporder=12
1885  case(181)
1886    brvsb="P"; intsb="6_4 2 2"; schsb="D6^5"; sporder=12
1887  case(182)
1888    brvsb="P"; intsb="6_3 2 2"; schsb="D6^6"; sporder=12
1889  case(183)
1890    brvsb="P"; intsb="6 m m"; schsb="C6v^1"; sporder=12
1891  case(184)
1892    brvsb="P"; intsb="6 c c"; schsb="C6v^2"; sporder=12
1893  case(185)
1894    brvsb="P"; intsb="6_3 c m"; schsb="C6v^3"; sporder=12
1895  case(186)
1896    brvsb="P"; intsb="6_3 m c"; schsb="C6v^4"; sporder=12
1897  case(187)
1898    brvsb="P"; intsb="-6 m 2"; schsb="D3h^1"; sporder=12
1899  case(188)
1900    brvsb="P"; intsb="-6 c 2"; schsb="D3h^2"; sporder=12
1901  case(189)
1902    brvsb="P"; intsb="-6 2 m"; schsb="D3h^3"; sporder=12
1903  case(190)
1904    brvsb="P"; intsb="-6 2 c"; schsb="D3h^4"; sporder=12
1905  case(191)
1906    brvsb="P"; intsb="6/m m m"; schsb="D6h^1"; sporder=24
1907  case(192)
1908    brvsb="P"; intsb="6/m c c"; schsb="D6h^2"; sporder=24
1909  case(193)
1910    brvsb="P"; intsb="6_3/m c m"; schsb="D6h^3"; sporder=24
1911  case(194)
1912    brvsb="P"; intsb="6_3/m m c"; schsb="D6h^4"; sporder=24
1913  case(195)
1914    brvsb="P"; intsb="2 3"; schsb="T^1"; sporder=12
1915  case(196)
1916    brvsb="F"; intsb="2 3"; schsb="T^2"; sporder=48
1917  case(197)
1918    brvsb="I"; intsb="2 3"; schsb="T^3"; sporder=24
1919  case(198)
1920    brvsb="P"; intsb="2_1 3"; schsb="T^4"; sporder=12
1921  case(199)
1922    brvsb="I"; intsb="2_1 3"; schsb="T^5"; sporder=24
1923  case(200)
1924    brvsb="P"; intsb="m -3"; schsb="Th^1"; sporder=24
1925  case(201)
1926    brvsb="P"; intsb="n -3"; schsb="Th^2"; sporder=24
1927    if (spgorig==1) then
1928      intsbl="n -3 _1"
1929    else if (spgorig==2)then
1930      intsbl="n -3 _2"
1931    else
1932      intsbl="intsbl to be determined"
1933    end if
1934  case(202)
1935    brvsb="F"; intsb="m -3"; schsb="Th^3"; sporder=96
1936  case(203)
1937    brvsb="F"; intsb="d -3"; schsb="Th^4"; sporder=96
1938    if (spgorig==1) then
1939      intsbl="d -3 _1"
1940    else if (spgorig==2) then
1941      intsbl="d -3 _2"
1942    else
1943      intsbl="intsbl to be determined"
1944    end if
1945  case(204)
1946    brvsb="I"; intsb="m -3"; schsb="Th^5"; sporder=48
1947  case(205)
1948    brvsb="P"; intsb="a -3"; schsb="Th^6"; sporder=24
1949  case(206)
1950    brvsb="I"; intsb="a -3"; schsb="Th^7"; sporder=48
1951  case(207)
1952    brvsb="P"; intsb="4 3 2"; schsb="O^1"; sporder=24
1953  case(208)
1954    brvsb="P"; intsb="4_2 3 2"; schsb="O^2"; sporder=24
1955  case(209)
1956    brvsb="F"; intsb="4 3 2"; schsb="O^3"; sporder=96
1957  case(210)
1958    brvsb="F"; intsb="4_1 3 2"; schsb="O^4"; sporder=96
1959  case(211)
1960    brvsb="I"; intsb="4 3 2"; schsb="O^5"; sporder=48
1961  case(212)
1962    brvsb="P"; intsb="4_3 3 2"; schsb="O^6"; sporder=24
1963  case(213)
1964    brvsb="P"; intsb="4_1 3 2"; schsb="O^7"; sporder=24
1965  case(214)
1966    brvsb="I"; intsb="4_1 3 2"; schsb="O^8"; sporder=48
1967  case(215)
1968    brvsb="P"; intsb="-4 3 m"; schsb="Td^1"; sporder=24
1969  case(216)
1970    brvsb="F"; intsb="-4 3 m"; schsb="Td^2"; sporder=96
1971  case(217)
1972    brvsb="I"; intsb="-4 3 m"; schsb="Td^3"; sporder=48
1973  case(218)
1974    brvsb="P"; intsb="-4 3 n"; schsb="Td^4"; sporder=24
1975  case(219)
1976    brvsb="F"; intsb="-4 3 c"; schsb="Td^5"; sporder=96
1977  case(220)
1978    brvsb="I"; intsb="-4 3 d"; schsb="Td^6"; sporder=48
1979  case(221)
1980    brvsb="P"; intsb="m -3 m"; schsb="Oh^1"; sporder=48
1981  case(222)
1982    brvsb="P"; intsb="n -3 n"; schsb="Oh^2"; sporder=48
1983    if (spgorig==1) then
1984      intsbl="n -3 n _1"
1985    else if (spgorig==2) then
1986      intsbl="n -3 n _2"
1987    else
1988      intsbl="intsbl to be determined"
1989    end if
1990  case(223)
1991    brvsb="P"; intsb="m -3 n"; schsb="Oh^3"; sporder=48
1992  case(224)
1993    brvsb="P"; intsb="n -3 m"; schsb="Oh^4"; sporder=48
1994    if (spgorig==1) then
1995      intsbl="n -3 m _1"
1996    else if (spgorig==2)then
1997      intsbl="n -3 m _2"
1998    else
1999      intsbl="intsbl to be determined"
2000    end if
2001  case(225)
2002    brvsb="F"; intsb="m -3 m"; schsb="Oh^5"; sporder=192
2003  case(226)
2004    brvsb="F"; intsb="m -3 c"; schsb="Oh^6"; sporder=192
2005  case(227)
2006    brvsb="F"; intsb="d -3 m"; schsb="Oh^7"; sporder=192
2007    if (spgorig==1) then
2008      intsbl="d -3 m _1"
2009    else if (spgorig==2) then
2010      intsbl="d -3 m _2"
2011    else
2012      intsbl="intsbl to be determined"
2013    end if
2014  case(228)
2015    brvsb="F"; intsb="d -3 c"; schsb="Oh^8"; sporder=192
2016    if (spgorig==1) then
2017      intsbl="d -3 c _1"
2018    else if (spgorig==2) then
2019      intsbl="d -3 c _2"
2020    else
2021      intsbl="intsbl to be determined"
2022    end if
2023  case(229)
2024    brvsb="I"; intsb="m -3 m"; schsb="Oh^9"; sporder=96
2025  case(230)
2026    brvsb="I"; intsb="a -3 d"; schsb="Oh^10"; sporder=96
2027  end select
2028 
2029  if(trim(intsbl)=="same")intsbl=intsb
2030 
2031 !Assignment of the point group number
2032  if(spgroup<=2)then  ! Triclinic system
2033    select case(spgroup)
2034    case (1)
2035      ptintsb="1"; ptschsb="C1"
2036    case (2)
2037      ptintsb="-1"; ptschsb="Ci"
2038    end select
2039  else if(spgroup<=15)then  ! Monoclinic system
2040    select case(spgroup)
2041    case (3:5)
2042      ptintsb="2"; ptschsb="C2"
2043    case (6:9)
2044      ptintsb="m"; ptschsb="Cs = C1h "
2045    case (10:15)
2046      ptintsb="2/m"; ptschsb="C2h"
2047    end select
2048  else if(spgroup<=74)then  ! Orthorhombic system
2049    select case(spgroup)
2050    case (16:24)
2051      ptintsb="2 2 2"; ptschsb="D2"
2052    case (25:46)
2053      ptintsb="m m 2"; ptschsb="C2v"
2054    case (47:74)
2055      ptintsb="m m m"; ptschsb="D2h"
2056    end select
2057  else if(spgroup<=142)then  ! Tetragonal system
2058    select case(spgroup)
2059    case (75:80)
2060      ptintsb="4"; ptschsb="C4"
2061    case (81,82)
2062      ptintsb="-4"; ptschsb="S4"
2063    case (83:88)
2064      ptintsb="4/m"; ptschsb="C4h"
2065    case (89:98)
2066      ptintsb="4 2 2"; ptschsb="D4"
2067    case (99:110)
2068      ptintsb="4 m m"; ptschsb="C4v"
2069    case (111:114,121,122)
2070      ptintsb="-4 2 m"; ptschsb="D2d^1"
2071    case (115:120)
2072      ptintsb="-4 m 2"; ptschsb="D2h^2"
2073    case (123:142)
2074      ptintsb="4/m m m"
2075      ptschsb="D4h"
2076    end select
2077  else if(spgroup<=167)then  ! Trigonal system
2078    select case(spgroup)
2079    case (143:146)
2080      ptintsb="3"; ptschsb="C3"
2081    case (147,148)
2082      ptintsb="-3"; ptschsb="C3i"
2083    case (149,151,153)
2084      ptintsb="3 1 2"; ptschsb="D3^1"
2085    case (150,152,154,155)
2086      ptintsb="3 2 1"; ptschsb="D3^2"
2087    case (156,158,160,161)
2088      ptintsb="3 m 1"; ptschsb="C3v^1"
2089    case (157,159)
2090      ptintsb="3 1 m"; ptschsb="C3v^2"
2091    case (162,163)
2092      ptintsb="-3 1 m"; ptschsb="D3d^1"
2093    case (164:167)
2094      ptintsb="-3 m 1"; ptschsb="D3d^2"
2095    end select
2096  else if(spgroup<=194)then  ! Hexagonal system
2097    select case(spgroup)
2098    case (168:173)
2099      ptintsb="6"; ptschsb="C6"
2100    case (174)
2101      ptintsb="-6"; ptschsb="C3h"
2102    case (175,176)
2103      ptintsb="6/m"; ptschsb="C6h"
2104    case (177:182)
2105      ptintsb="6 2 2"; ptschsb="D6"
2106    case (183:186)
2107      ptintsb="6 m m"; ptschsb="C6v"
2108    case (187,188)
2109      ptintsb="-6 m 2"; ptschsb="D3h^1"
2110    case (189,190)
2111      ptintsb="-6 2 m"; ptschsb="D3h^2"
2112    case (191:194)
2113      ptintsb="6/m m m"; ptschsb="D6h"
2114    end select
2115  else                        ! Cubic system
2116    select case(spgroup)
2117    case (195:199)
2118      ptintsb="2 3"; ptschsb="T"
2119    case (200:206)
2120      ptintsb="m 3"; ptschsb="Th"
2121    case (207:214)
2122      ptintsb="4 3 2"; ptschsb="O"
2123    case (215:220)
2124      ptintsb="4 3 m"; ptschsb="Td"
2125    case (221:230)
2126      ptintsb="m -3 m"; ptschsb="Oh"
2127    end select
2128  end if
2129 
2130 end subroutine spgdata

m_spgdata/symptgroup [ Functions ]

[ Top ] [ m_spgdata ] [ Functions ]

NAME

 symptgroup

FUNCTION

 Derive the name of the point group (+holohedry), from symrel.
 Warning: might have to change the holohedry hR to hP, if hexagonal axes

INPUTS

 nsym=actual number of symmetries
 symrel(3,3,nsym)=nsym symmetry operations in real space in terms of primitive translations

OUTPUT

 iholohedry=holohedry number
 ptgroup=symmetry point group

PARENTS

      symanal,symbrav

CHILDREN

      symdet

SOURCE

2476 subroutine symptgroup(iholohedry,nsym,ptgroup,symrel)
2477 
2478 
2479 !This section has been created automatically by the script Abilint (TD).
2480 !Do not modify the following lines by hand.
2481 #undef ABI_FUNC
2482 #define ABI_FUNC 'symptgroup'
2483 !End of the abilint section
2484 
2485  implicit none
2486 
2487 !Arguments ------------------------------------
2488 !scalars
2489  integer,intent(in) :: nsym
2490  integer,intent(out) :: iholohedry
2491  character(len=5),intent(out) :: ptgroup
2492 !arrays
2493  integer,intent(in) :: symrel(3,3,nsym)
2494 
2495 !Local variables-------------------------------
2496 !scalars
2497  integer :: inversion,iorder,isym
2498  character(len=500) :: message
2499 !arrays
2500  integer :: identity(3,3),matrix(3,3),n_axes(-6:6),trial(3,3)
2501  integer,allocatable :: determinant(:),order(:),root_invers(:)
2502  character(len=2),allocatable :: ptsym(:)
2503 
2504 !**************************************************************************
2505 
2506 !DEBUG
2507 !write(std_out,*)' symptgroup : enter'
2508 !do isym=1,nsym
2509 !write(std_out,'(i3,2x,9i3)' )isym,symrel(:,:,isym)
2510 !end do
2511 !ENDDEBUG
2512 
2513  identity(:,:)=0
2514  identity(1,1)=1 ; identity(2,2)=1 ; identity(3,3)=1
2515  n_axes(:)=0
2516 
2517  ABI_ALLOCATE(determinant,(nsym))
2518  ABI_ALLOCATE(order,(nsym))
2519  ABI_ALLOCATE(ptsym,(nsym))
2520  ABI_ALLOCATE(root_invers,(nsym))
2521 
2522 !Get the determinant
2523  call symdet(determinant,nsym,symrel)
2524 
2525 !Get the order of each the symmetry operation, as well as the maximal order
2526 !Also, examine whether each symmetry operation is the inversion, or a root
2527 !of the inversion (like -3)
2528 !Finally, decide which kind of point symmetry operation it is
2529  do isym=1,nsym
2530 
2531    trial(:,:)=identity(:,:)
2532    matrix(:,:)=symrel(:,:,isym)
2533    order(isym)=0
2534    root_invers(isym)=0
2535    do iorder=1,6
2536      trial=matmul(matrix,trial)
2537      if(sum((trial-identity)**2)==0)then
2538        order(isym)=iorder
2539        exit
2540      end if
2541      if(sum((trial+identity)**2)==0)then
2542        root_invers(isym)=iorder
2543        if(iorder==1)inversion=isym
2544      end if
2545    end do
2546    if(order(isym)==0)then
2547      write(message, '(a,i0,a)' )' The symmetry operation number',isym,' is not a root of unity'
2548      MSG_BUG(message)
2549    end if
2550 
2551 !  determinant, order and root_invers are enough to determine the
2552 !  kind of symmetry operation
2553    ptsym(isym)='no'
2554    select case(order(isym))
2555    case(1)
2556      ptsym(isym)=' 1' ; n_axes(1)=n_axes(1)+1
2557    case(2)
2558      if(determinant(isym)== 1)then
2559        ptsym(isym)=' 2' ; n_axes(2)=n_axes(2)+1
2560      else if(determinant(isym)==-1 .and. root_invers(isym)==1)then
2561        ptsym(isym)='-1' ; n_axes(-1)=n_axes(-1)+1
2562      else if(determinant(isym)==-1 .and. root_invers(isym)==0)then
2563        ptsym(isym)='-2' ; n_axes(-2)=n_axes(-2)+1
2564      end if
2565    case(3)
2566      ptsym(isym)=' 3' ; n_axes(3)=n_axes(3)+1
2567    case(4)
2568      if(determinant(isym)== 1)then
2569        ptsym(isym)=' 4' ; n_axes(4)=n_axes(4)+1
2570      else if(determinant(isym)==-1)then
2571        ptsym(isym)='-4' ; n_axes(-4)=n_axes(-4)+1
2572      end if
2573    case(6)
2574      if(determinant(isym)== 1)then
2575        ptsym(isym)=' 6' ; n_axes(6)=n_axes(6)+1
2576      else if(determinant(isym)==-1 .and. root_invers(isym)==3)then
2577        ptsym(isym)='-3' ; n_axes(-3)=n_axes(-3)+1
2578      else if(determinant(isym)==-1 .and. root_invers(isym)==0)then
2579        ptsym(isym)='-6' ; n_axes(-6)=n_axes(-6)+1
2580      end if
2581    end select
2582 
2583    if(ptsym(isym)=='no')then
2584      write(message,'(a,i4,a,a,a,i4,a,a,i4,a,a,i4)' )&
2585 &     'The symmetry operation number',isym,' could not be identified',ch10,&
2586 &     'order(isym)      =',order(isym),ch10,&
2587 &     'determinant(isym)=',determinant(isym),ch10,&
2588 &     'root_invers(isym)=',root_invers(isym)
2589      MSG_BUG(message)
2590    end if
2591 
2592  end do
2593 
2594  iholohedry=0
2595  if     (sum((n_axes-(/0,0,0,0,0,0, 0 ,1,0,0,0,0,0/))**2)==0)then
2596    ptgroup='    1' ; iholohedry=1
2597  else if(sum((n_axes-(/0,0,0,0,0,1, 0 ,1,0,0,0,0,0/))**2)==0)then
2598    ptgroup='   -1' ; iholohedry=1
2599 
2600  else if(sum((n_axes-(/0,0,0,0,0,0, 0 ,1,1,0,0,0,0/))**2)==0)then
2601    ptgroup='    2' ; iholohedry=2
2602  else if(sum((n_axes-(/0,0,0,0,1,0, 0 ,1,0,0,0,0,0/))**2)==0)then
2603    ptgroup='   -2' ; iholohedry=2
2604  else if(sum((n_axes-(/0,0,0,0,1,1, 0 ,1,1,0,0,0,0/))**2)==0)then
2605    ptgroup='  2/m' ; iholohedry=2
2606 
2607  else if(sum((n_axes-(/0,0,0,0,0,0, 0 ,1,3,0,0,0,0/))**2)==0)then
2608    ptgroup='  222' ; iholohedry=3
2609  else if(sum((n_axes-(/0,0,0,0,2,0, 0 ,1,1,0,0,0,0/))**2)==0)then
2610    ptgroup='  mm2' ; iholohedry=3
2611  else if(sum((n_axes-(/0,0,0,0,3,1, 0 ,1,3,0,0,0,0/))**2)==0)then
2612    ptgroup='  mmm' ; iholohedry=3
2613 
2614  else if(sum((n_axes-(/0,0,0,0,0,0, 0 ,1,1,0,2,0,0/))**2)==0)then
2615    ptgroup='    4' ; iholohedry=4
2616  else if(sum((n_axes-(/0,0,2,0,0,0, 0 ,1,1,0,0,0,0/))**2)==0)then
2617    ptgroup='   -4' ; iholohedry=4
2618  else if(sum((n_axes-(/0,0,2,0,1,1, 0 ,1,1,0,2,0,0/))**2)==0)then
2619    ptgroup='  4/m' ; iholohedry=4
2620  else if(sum((n_axes-(/0,0,0,0,0,0, 0 ,1,5,0,2,0,0/))**2)==0)then
2621    ptgroup='  422' ; iholohedry=4
2622  else if(sum((n_axes-(/0,0,0,0,4,0, 0 ,1,1,0,2,0,0/))**2)==0)then
2623    ptgroup='  4mm' ; iholohedry=4
2624  else if(sum((n_axes-(/0,0,2,0,2,0, 0 ,1,3,0,0,0,0/))**2)==0)then
2625    ptgroup=' -42m' ; iholohedry=4
2626  else if(sum((n_axes-(/0,0,2,0,5,1, 0 ,1,5,0,2,0,0/))**2)==0)then
2627    ptgroup='4/mmm' ; iholohedry=4
2628 
2629  else if(sum((n_axes-(/0,0,0,0,0,0, 0 ,1,0,2,0,0,0/))**2)==0)then
2630    ptgroup='    3' ; iholohedry=5
2631  else if(sum((n_axes-(/0,0,0,2,0,1, 0 ,1,0,2,0,0,0/))**2)==0)then
2632    ptgroup='   -3' ; iholohedry=5
2633  else if(sum((n_axes-(/0,0,0,0,0,0, 0 ,1,3,2,0,0,0/))**2)==0)then
2634    ptgroup='   32' ; iholohedry=5
2635  else if(sum((n_axes-(/0,0,0,0,3,0, 0 ,1,0,2,0,0,0/))**2)==0)then
2636    ptgroup='   3m' ; iholohedry=5
2637  else if(sum((n_axes-(/0,0,0,2,3,1, 0 ,1,3,2,0,0,0/))**2)==0)then
2638    ptgroup='  -3m' ; iholohedry=5
2639 
2640  else if(sum((n_axes-(/0,0,0,0,0,0, 0 ,1,1,2,0,0,2/))**2)==0)then
2641    ptgroup='    6' ; iholohedry=6
2642  else if(sum((n_axes-(/2,0,0,0,1,0, 0 ,1,0,2,0,0,0/))**2)==0)then
2643    ptgroup='   -6' ; iholohedry=6
2644  else if(sum((n_axes-(/2,0,0,2,1,1, 0 ,1,1,2,0,0,2/))**2)==0)then
2645    ptgroup='  6/m' ; iholohedry=6
2646  else if(sum((n_axes-(/0,0,0,0,0,0, 0 ,1,7,2,0,0,2/))**2)==0)then
2647    ptgroup='  622' ; iholohedry=6
2648  else if(sum((n_axes-(/0,0,0,0,6,0, 0 ,1,1,2,0,0,2/))**2)==0)then
2649    ptgroup='  6mm' ; iholohedry=6
2650  else if(sum((n_axes-(/2,0,0,0,4,0, 0 ,1,3,2,0,0,0/))**2)==0)then
2651    ptgroup=' -62m' ; iholohedry=6
2652  else if(sum((n_axes-(/2,0,0,2,7,1, 0 ,1,7,2,0,0,2/))**2)==0)then
2653    ptgroup='6/mmm' ; iholohedry=6
2654 
2655  else if(sum((n_axes-(/0,0,0,0,0,0, 0 ,1,3,8,0,0,0/))**2)==0)then
2656    ptgroup='   23' ; iholohedry=7
2657  else if(sum((n_axes-(/0,0,0,8,3,1, 0 ,1,3,8,0,0,0/))**2)==0)then
2658    ptgroup='  m-3' ; iholohedry=7
2659  else if(sum((n_axes-(/0,0,0,0,0,0, 0 ,1,9,8,6,0,0/))**2)==0)then
2660    ptgroup='  432' ; iholohedry=7
2661  else if(sum((n_axes-(/0,0,6,0,6,0, 0 ,1,3,8,0,0,0/))**2)==0)then
2662    ptgroup=' -43m' ; iholohedry=7
2663  else if(sum((n_axes-(/0,0,6,8,9,1, 0 ,1,9,8,6,0,0/))**2)==0)then
2664    ptgroup=' m-3m' ; iholohedry=7
2665 
2666  end if
2667 
2668  if(iholohedry==0)then
2669    MSG_ERROR_CLASS('Could not find the point group', "TolSymError")
2670  end if
2671 
2672 !DEBUG
2673 !do isym=1,nsym
2674 !write(std_out,'(a,3i5)' )&
2675 !&  ' symptgroup : isym,determinant,order=',isym,determinant(isym),order(isym)
2676 !end do
2677 
2678 !write(std_out,'(a,13i3)' )' symptgroup : n_axes(-6:6)=',n_axes(-6:6)
2679 !write(std_out,*)' iholohedry, ptgroup=',iholohedry,',',ptgroup
2680 !ENDDEBUG
2681 
2682  ABI_DEALLOCATE(determinant)
2683  ABI_DEALLOCATE(order)
2684  ABI_DEALLOCATE(ptsym)
2685  ABI_DEALLOCATE(root_invers)
2686 
2687 end subroutine symptgroup