TABLE OF CONTENTS


ABINIT/m_symsg [ Modules ]

[ Top ] [ Modules ]

NAME

  m_symsg

FUNCTION

COPYRIGHT

  Copyright (C) 1999-2018 ABINIT group (RC,XG)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

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_symsg
28 
29  use defs_basis
30  use m_errors
31  use m_abicore
32 
33  use m_time,     only : timab
34  use m_spgdata,  only : spgdata
35 
36  implicit none
37 
38  private

m_symsg/bldgrp [ Functions ]

[ Top ] [ m_symsg ] [ Functions ]

NAME

 bldgrp

FUNCTION

 Yields all the symmetry operations starting from the generators.
 Applies all the generators onto themselves, and obtains all the other operations.
 Iterates until it reaches nsym.

INPUTS

 msym = default number of symmetry operations
 nsym = number of symmetry operations
 symafm(msym)=(anti)ferromagnetic part of symmetry operations
 symrel(3,3,msym) = 3D matrix containg symmetry operations
 tnons(3,msym) = 2D matrix containing translations of the symmery operations

OUTPUT

 symafm(msym)=(anti)ferromagnetic part of symmetry operations
 symrel(3,3,msym) = 3D matrix containg symmetry operations
 tnons(3,msym) = 2D matrix containing translations of the symmery operations

SIDE EFFECTS

 nogen = number of generators, number of operations to be applied onto themselves

PARENTS

      symsgcube,symsghexa,symsgortho,symsgtetra

CHILDREN

SOURCE

2436 subroutine bldgrp(msym,nogen,nsym,symafm,symrel,tnons)
2437 
2438 
2439 !This section has been created automatically by the script Abilint (TD).
2440 !Do not modify the following lines by hand.
2441 #undef ABI_FUNC
2442 #define ABI_FUNC 'bldgrp'
2443 !End of the abilint section
2444 
2445  implicit none
2446 
2447 !Arguments ------------------------------------
2448 !scalars
2449  integer,intent(in) :: msym,nsym
2450  integer,intent(inout) :: nogen
2451 !arrays
2452  integer,intent(inout) :: symafm(msym),symrel(3,3,msym)
2453  real(dp),intent(inout) :: tnons(3,msym)
2454 
2455 !Local variables ------------------------------
2456 !matrintoper(3,3) & matrinttransl(3) are intermediate arrays of the new
2457 !      symmetry operations obtained, in order to check their uniqueness.
2458 !flagop,flagtr = flags used during the checking of the similarity between
2459 !      the obtained operation and the already existent ones
2460 !ii,ijk,ijkl,jjj,kk = counters in the cycles
2461 !scalars
2462  integer :: flagop,flagtr,ii,ijk,ijkl,jjj,kk,matrintsymafm,nogen_new
2463  real(dp) :: nastyzero
2464  character(len=500) :: message
2465 !arrays
2466  integer :: bcksymafm(2*msym),bcksymrel(3,3,2*msym),matrintoper(3,3)
2467  real(dp) :: bcktnons(3,2*msym),matrinttransl(3)
2468 
2469 ! *************************************************************************
2470 
2471  nastyzero=0.1
2472 
2473 !DEBUG
2474 ! write(std_out,*)' bldgrp : enter, builds the space group symmetry '
2475 ! write(std_out,*)' bldgrp : number of generators : ',nogen
2476 ! write(std_out,*)' bldgrp : nsym,msym=',nsym,msym
2477 !ENDDEBUG
2478 
2479  if (nogen<1) then
2480    write(message, '(a,i4,a,a,a,a,a)' )&
2481 &   'The number of generators nogen is ',nogen,&
2482 &   'and it should be greater than one',ch10,&
2483 &   'This is not allowed.  ',ch10,&
2484 &   'Action: Contact ABINIT group '
2485    MSG_ERROR(message)
2486  end if
2487 
2488 !Transfer the generators to bcksymrel
2489  do ii=1,nogen
2490    bcksymrel(:,:,ii)=symrel(:,:,ii)
2491    bcktnons(:,ii)=tnons(:,ii)
2492    bcksymafm(ii)=symafm(ii)
2493  end do
2494 
2495 !DEBUG
2496  write(std_out,*)' Describe the different generators (index,symrel,tnons,symafm)'
2497  do ii=1,nogen
2498    write(std_out,'(i3,2x,9i3,3es12.2,i3)')ii,symrel(:,:,ii),tnons(:,ii),symafm(ii)
2499  end do
2500 !ENDDEBUG
2501 
2502 !Simply iterate until the group is complete
2503  do ijkl=1,nsym
2504 
2505 !  DEBUG
2506 !  write(std_out,*)' bldgrp : in loop, ijkl,nogen=',ijkl,nogen
2507 !  ENDDEBUG
2508 
2509    nogen_new=nogen
2510 
2511    do jjj=2,nogen
2512      do kk=2,nogen
2513 
2514 !      Computing block of the new symmetry operation according to:
2515 !      !   $ { R1 | v1 }{ R2 | v2 } = { R1.R2 | v1+R1.v2 } $
2516        matrintoper(:,:) = matmul(bcksymrel(:,:,jjj),bcksymrel(:,:,kk))
2517        matrinttransl(:) = bcktnons(:,jjj)+matmul(bcksymrel(:,:,jjj),bcktnons(:,kk))
2518        matrintsymafm    = bcksymafm(jjj)*bcksymafm(kk)
2519 
2520 !      Rescaling translation between 0 and 1
2521        do ii=1,3
2522          if (matrinttransl(ii)>=0.9) then
2523            do while (matrinttransl(ii)>=0.9)
2524              matrinttransl(ii)=matrinttransl(ii)-1.0
2525            end do
2526          end if
2527          if (matrinttransl(ii)<0.0) then
2528            do while (matrinttransl(ii)<0.0)
2529              matrinttransl(ii)=matrinttransl(ii)+1.0
2530            end do
2531          end if
2532          if ( abs(matrinttransl(ii))<nastyzero) matrinttransl(ii)=0.0
2533          if ( abs(matrinttransl(ii)-1.0)<nastyzero) matrinttransl(ii)=0.0
2534        end do
2535 
2536 !      Cheking block to validate the new symmetry operation
2537        do ijk=1,nogen_new
2538 
2539          flagop=0 ; flagtr=0
2540 
2541 !        Check for rotation similarity
2542          if(sum((matrintoper-bcksymrel(:,:,ijk))**2)==0)flagop=1
2543 
2544 !        Check for translation similarity
2545          if(maxval((matrinttransl-bcktnons(:,ijk))**2)<nastyzero**2)flagtr=1
2546 
2547          if(flagop+flagtr==2)exit
2548 
2549        end do
2550 
2551 !      Add the new determined symmetry if it is unique
2552        if (flagtr+flagop<2) then
2553          nogen_new=nogen_new+1
2554 !        DEBUG
2555 !         write(std_out,*)' added one more symmetry : nogen_new=',nogen_new
2556 !        ENDDEBUG
2557          bcksymrel(:,:,nogen_new)=matrintoper(:,:)
2558          bcktnons(:,nogen_new)=matrinttransl(:)
2559          bcksymafm(nogen_new)=matrintsymafm
2560        end if
2561 
2562      end do
2563    end do
2564 
2565    nogen=nogen_new
2566 
2567    if(nogen==nsym)exit
2568 
2569  end do
2570 
2571 !Transfer of the calculated symmetry to the routine output
2572  if (nogen==nsym) then
2573    symrel(:,:,1:nsym)=bcksymrel(:,:,1:nsym)
2574    tnons(:,1:nsym)=bcktnons(:,1:nsym)
2575    symafm(1:nsym)=bcksymafm(1:nsym)
2576  else
2577 !  Problem with the generation of the symmetry operations
2578    write(message, '(a,i7,a,a,i7)' )&
2579 &   'The symmetries obtained are  ',nogen,ch10,&
2580 &   'and they should be ',nsym
2581    MSG_BUG(message)
2582  end if
2583 
2584 !DEBUG
2585 !write(std_out,*)' bldgrp : exit with  ',nogen,' operation symmetries'
2586 !ENDDEBUG
2587 
2588 end subroutine bldgrp

m_symsg/symsgcube [ Functions ]

[ Top ] [ m_symsg ] [ Functions ]

NAME

 symsgcube

FUNCTION

 Generate all the symmetry operations starting from the space group symbol
 for the cubic groups (according to the International Tables of Crystallography, 1983)

INPUTS

 msym = default number of symmetries
 shubnikov= magnetic type of the space group to be generated
 spgaxor = the possible orientation of the axes system
 spgorig = the origin choice (1 or 2) for the axes system
 spgroup = the numeric symbol of the space groups
 spgroupma = number of the magnetic space group

OUTPUT

 nsym = the number of symmetry operations
 symafm(nsym)=(anti)ferromagnetic part of symmetry operations
 symrel(3,3,nsym) = 3D matrix containg symmetry operations
 tnons(3,nsym) = 2D matrix containing translations associated

PARENTS

      gensymspgr

CHILDREN

      bldgrp,spgdata,timab

SOURCE

 81 subroutine symsgcube(msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma,symafm,symrel,tnons)
 82 
 83 
 84 !This section has been created automatically by the script Abilint (TD).
 85 !Do not modify the following lines by hand.
 86 #undef ABI_FUNC
 87 #define ABI_FUNC 'symsgcube'
 88 !End of the abilint section
 89 
 90  implicit none
 91 
 92 !Arguments ------------------------------------
 93 !scalars
 94  integer,intent(in) :: msym,shubnikov,spgaxor,spgorig,spgroup,spgroupma
 95  integer,intent(inout) :: nsym !vz_i
 96 !arrays
 97  integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i
 98  real(dp),intent(out) :: tnons(3,msym)
 99 
100 !Local variables -----------------------------
101 !scalars
102  integer :: ii,nogen,sporder
103  character(len=1) :: brvsb
104  character(len=15) :: intsb,ptintsb,ptschsb,schsb
105  character(len=35) :: intsbl
106 !arrays
107  integer :: gen1(3,3),gen2(3,3),gen3(3,3),gen4(3,3),gen5(3,3),gen6(3,3)
108  integer :: gen7(3,3),gen8(3,3),gen9(3,3),genmmm(3,3),genmmp(3,3),genmpm(3,3)
109  integer :: genmpp(3,3),genpmm(3,3),genpmp(3,3),genppm(3,3),genrot(3,3)
110  integer :: genswm(3,3),genswp(3,3)
111  real(dp) :: tsec(2)
112 
113 !*************************************************************************
114 
115 !DEBUG
116 !write(std_out,*) ' symsgcube : enter with spgroup ',spgroup,' and ',' origin choice ',spgorig
117 !write(std_out,*) ' msym,nsym=',msym,nsym
118 !ENDDEBUG
119 
120 
121 !The identity operation belongs to all space groups
122  symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1
123 
124 !Initialize the associated translations matrix to 0
125  do ii=1,msym
126    tnons(:,ii)= 0.0d0
127  end do
128  nogen=0
129 
130 !Predefine some generators
131  genmpp(:,:)=0 ; genmpp(1,1)=-1 ; genmpp(2,2)= 1 ; genmpp(3,3)= 1
132  genpmp(:,:)=0 ; genpmp(1,1)= 1 ; genpmp(2,2)=-1 ; genpmp(3,3)= 1
133  genppm(:,:)=0 ; genppm(1,1)= 1 ; genppm(2,2)= 1 ; genppm(3,3)=-1
134  genpmm(:,:)=0 ; genpmm(1,1)= 1 ; genpmm(2,2)=-1 ; genpmm(3,3)=-1
135  genmpm(:,:)=0 ; genmpm(1,1)=-1 ; genmpm(2,2)= 1 ; genmpm(3,3)=-1
136  genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)= 1
137  genmmm(:,:)=0 ; genmmm(1,1)=-1 ; genmmm(2,2)=-1 ; genmmm(3,3)=-1
138  genrot(:,:)=0 ; genrot(1,3)=1 ; genrot(3,2)=1 ; genrot(2,1)=1  !reshape((/0,0,1,1,0,0,0,1,0/),(/3,3/),(/0,0/),(/2,1/) )
139  genswm(:,:)=0 ; genswm(2,1)=1 ; genswm(1,2)=1 ; genswm(3,3)=-1 !reshape((/0,1,0,1,0,0,0,0,-1/),(/3,3/),(/0,0/),(/2,1/) )
140  genswp(:,:)=0 ; genswp(2,1)=1 ; genswp(1,2)=1 ; genswp(3,3)=1  !reshape((/0,1,0,1,0,0,0,0,1/),(/3,3/),(/0,0/),(/2,1/) )
141 
142 !Because of problems with the IBM compiler, that does not like reshape
143 !operations, define 9 basic matrices
144  gen1(:,:)=0 ; gen1(1,1)=1
145  gen2(:,:)=0 ; gen2(1,2)=1
146  gen3(:,:)=0 ; gen3(1,3)=1
147  gen4(:,:)=0 ; gen4(2,1)=1
148  gen5(:,:)=0 ; gen5(2,2)=1
149  gen6(:,:)=0 ; gen6(2,3)=1
150  gen7(:,:)=0 ; gen7(3,1)=1
151  gen8(:,:)=0 ; gen8(3,2)=1
152  gen9(:,:)=0 ; gen9(3,3)=1
153 
154 !Default non-magnetic behaviour
155  symafm(1:msym)=1
156 
157 !*************************************************************************
158 
159 !Treat CUBIC groups
160 
161  if(195<=spgroup .and. spgroup<=230)then
162 
163    select case(spgroup)
164    case (195,196,197)        !P23, F23, I23
165      symrel(:,:,2) = genrot(:,:)
166      symrel(:,:,3) = genmmp(:,:)
167      symrel(:,:,4) = genmpm(:,:)
168      symrel(:,:,5) = genpmm(:,:)
169      nogen=5
170    case (200,202,204)                !PmB3, FmB3, ImB3
171      symrel(:,:,2) = genrot(:,:)
172      symrel(:,:,3) = genmmp(:,:)
173      symrel(:,:,4) = genmpm(:,:)
174      nogen=4
175    case (198,199)                !P213, I213
176      symrel(:,:,2) = genrot(:,:)
177      symrel(:,:,3) = genmmp(:,:)
178      tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
179      symrel(:,:,4) = genmpm(:,:)
180      tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
181      symrel(:,:,5) = genpmm(:,:)
182      tnons(:,5)=(/0.5d0,0.5d0,0.d0/)
183      symrel(:,:,6) = gen3(:,:)-gen4(:,:)-gen8(:,:)
184      tnons(:,6)=(/0.5d0,0.5d0,0.d0/)
185      symrel(:,:,7) = -gen3(:,:)-gen4(:,:)+gen8(:,:)
186      tnons(:,7)=(/0.5d0,0.d0,0.5d0/)
187      symrel(:,:,8) = -gen3(:,:)+gen4(:,:)-gen8(:,:)
188      tnons(:,8)=(/0.d0,0.5d0,0.5d0/)
189      symrel(:,:,10) = -gen2(:,:)+gen6(:,:)-gen7(:,:)
190      tnons(:,10)=(/0.d0,0.5d0,0.5d0/)
191      symrel(:,:,11) =  gen2(:,:)-gen6(:,:)-gen7(:,:)
192      tnons(:,11)=(/0.5d0,0.5d0,0.0d0/)
193      symrel(:,:,12) = -gen2(:,:)-gen6(:,:)+gen7(:,:)
194      tnons(:,12)=(/0.5d0,0.d0,0.5d0/)
195      symrel(:,:,9) =  gen2(:,:)+gen6(:,:)+gen7(:,:)
196    case (201)                !Pn-3
197      if (spgorig==1) then
198        symrel(:,:,2) = genrot(:,:)
199        symrel(:,:,3) = genmmp(:,:)
200        symrel(:,:,4) = genmpm(:,:)
201        symrel(:,:,5) = genmmm(:,:)
202        tnons(:,5)=(/0.5d0,0.5d0,0.5d0/)
203        nogen=5
204        if(shubnikov==3)symafm(5)=-1
205        call bldgrp(msym,nogen,nsym,symafm,symrel,tnons)
206      else
207        symrel(:,:,2) = genrot(:,:)
208        symrel(:,:,3) = genmmp(:,:)
209        tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
210        symrel(:,:,4) = genmpm(:,:)
211        tnons(:,4)=(/0.5d0,0.d0,0.5d0/)
212        symrel(:,:,5) = genpmm(:,:)
213        tnons(:,5)=(/0.d0,0.5d0,0.5d0/)
214        symrel(:,:,6) =  gen3(:,:)-gen4(:,:)-gen8(:,:)
215        tnons(:,6)=(/0.d0,0.5d0,0.5d0/)
216        symrel(:,:,7) = -gen3(:,:)-gen4(:,:)+gen8(:,:)
217        tnons(:,7)=(/0.5d0,0.5d0,0.d0/)
218        symrel(:,:,8) = -gen3(:,:)+gen4(:,:)-gen8(:,:)
219        tnons(:,8)=(/0.5d0,0.d0,0.5d0/)
220        symrel(:,:,10) = -gen2(:,:)+gen6(:,:)-gen7(:,:)
221        tnons(:,10)=(/0.5d0,0.d0,0.5d0/)
222        symrel(:,:,11) =  gen2(:,:)-gen6(:,:)-gen7(:,:)
223        tnons(:,11)=(/0.d0,0.5d0,0.5d0/)
224        symrel(:,:,12) = -gen2(:,:)-gen6(:,:)+gen7(:,:)
225        tnons(:,12)=(/0.5d0,0.5d0,0.d0/)
226        symrel(:,:,9) =  gen2(:,:)+gen6(:,:)+gen7(:,:)
227        do ii=1,12
228          symrel(:,:,ii+12)= - symrel(:,:,ii)
229          tnons(:,ii+12)=tnons(:,ii)
230          if(shubnikov==3)symafm(ii+12)=-1
231        end do
232      end if
233      nogen=0
234    case (205,206)                !Pa-3,Ia-3
235      symrel(:,:,2) = genrot(:,:)
236      symrel(:,:,3) = genmmp(:,:)
237      tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
238      symrel(:,:,4) = genmpm(:,:)
239      tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
240      symrel(:,:,5) = genpmm(:,:)
241      tnons(:,5)=(/0.5d0,0.5d0,0.d0/)
242      symrel(:,:,6) =  gen3(:,:)-gen4(:,:)-gen8(:,:)
243      tnons(:,6)=(/0.5d0,0.5d0,0.d0/)
244      symrel(:,:,7) = -gen3(:,:)-gen4(:,:)+gen8(:,:)
245      tnons(:,7)=(/0.5d0,0.d0,0.5d0/)
246      symrel(:,:,8) = -gen3(:,:)+gen4(:,:)-gen8(:,:)
247      tnons(:,8)=(/0.0d0,0.5d0,0.5d0/)
248      symrel(:,:,9) =  gen2(:,:)+gen6(:,:)+gen7(:,:)
249      symrel(:,:,10) = -gen2(:,:)+gen6(:,:)-gen7(:,:)
250      tnons(:,10)=(/0.d0,0.5d0,0.5d0/)
251      symrel(:,:,11) =  gen2(:,:)-gen6(:,:)-gen7(:,:)
252      tnons(:,11)=(/0.5d0,0.5d0,0.d0/)
253      symrel(:,:,12) = -gen2(:,:)-gen6(:,:)+gen7(:,:)
254      tnons(:,12)=(/0.5d0,0.d0,0.5d0/)
255      nogen=0
256    case (203)                !FdB3
257      if (spgorig==1) then
258        symrel(:,:,2) = genrot(:,:)
259        symrel(:,:,3) = genmmp(:,:)
260        symrel(:,:,4) = genmpm(:,:)
261        nogen=4
262        nsym=12
263        call bldgrp(msym,nogen,nsym,symafm,symrel,tnons)
264        do ii=1,12
265          symrel(:,:,ii+12) = - symrel(:,:,ii)
266          tnons(:,ii+12)=(/0.25d0,0.25d0,0.25d0/)
267          if(shubnikov==3)symafm(ii+12)=-1
268        end do
269        nsym=24
270        nogen=0
271      else
272        symrel(:,:,2) = genrot(:,:)
273        symrel(:,:,3) = genmmp(:,:)
274        tnons(:,3)=(/0.25d0,0.25d0,0.d0/)
275        symrel(:,:,4) = genmpm(:,:)
276        tnons(:,4)=(/0.25d0,0.d0,0.25d0/)
277        symrel(:,:,5) = genpmm(:,:)
278        tnons(:,5)=(/0.d0,0.25d0,0.25d0/)
279        symrel(:,:,6) =  gen3(:,:)-gen4(:,:)-gen8(:,:)
280        tnons(:,6)=(/0.d0,0.25d0,0.25d0/)
281        symrel(:,:,7) = -gen3(:,:)-gen4(:,:)+gen8(:,:)
282        tnons(:,7)=(/0.25d0,0.25d0,0.d0/)
283        symrel(:,:,8) = -gen3(:,:)+gen4(:,:)-gen8(:,:)
284        tnons(:,8)=(/0.25d0,0.d0,0.25d0/)
285        symrel(:,:,10) = -gen2(:,:)+gen6(:,:)-gen7(:,:)
286        tnons(:,10)=(/0.25d0,0.d0,0.25d0/)
287        symrel(:,:,11) =  gen2(:,:)-gen6(:,:)-gen7(:,:)
288        tnons(:,11)=(/0.d0,0.25d0,0.25d0/)
289        symrel(:,:,12) = -gen2(:,:)-gen6(:,:)+gen7(:,:)
290        tnons(:,12)=(/0.25d0,0.25d0,0.d0/)
291        symrel(:,:,9) =  gen2(:,:)+gen6(:,:)+gen7(:,:)
292        do ii=1,12
293          symrel(:,:,ii+12) = - symrel(:,:,ii)
294          tnons(:,ii+12)=tnons(:,ii)
295          if(shubnikov==3)symafm(ii+12)=-1
296        end do
297      end if
298      nsym=24
299      nogen=0
300    case (207,209,211)        !P432, F432, I432, PmB3m, FmB3m, ImB3m
301      symrel(:,:,2) = genrot(:,:)
302      symrel(:,:,3) = genmmp(:,:)
303      symrel(:,:,4) = genmpm(:,:)
304      symrel(:,:,5) = genswm(:,:)
305      if(shubnikov==3)symafm(5)=-1
306      nogen=5
307    case (208)                !P4232
308      symrel(:,:,2) = genrot(:,:)
309      symrel(:,:,3) = genmmp(:,:)
310      symrel(:,:,4) = genmpm(:,:)
311      symrel(:,:,5) = genswm(:,:)
312      tnons(:,5)=(/0.5d0,0.5d0,0.5d0/)
313      if(shubnikov==3)symafm(5)=-1
314      nogen=5
315    case (210)                !F4132
316      symrel(:,:,2) = genrot(:,:)
317      symrel(:,:,3) = genmmp(:,:)
318      tnons(:,3)=(/0.d0,0.5d0,0.5d0/)
319      symrel(:,:,4) = genmpm(:,:)
320      tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
321      symrel(:,:,5) = genswm(:,:)
322      tnons(:,5)=(/0.75d0,0.25d0,0.75d0/)
323      if(shubnikov==3)symafm(5)=-1
324      nogen=5
325    case (212)                !P4332
326      symrel(:,:,2) = genrot(:,:)
327      symrel(:,:,3) = genmmp(:,:)
328      tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
329      symrel(:,:,4) = genmpm(:,:)
330      tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
331      symrel(:,:,5) = genswm(:,:)
332      tnons(:,5)=(/0.25d0,0.75d0,0.75d0/)
333      if(shubnikov==3)symafm(5)=-1
334      nogen=5
335    case (213,214)                !P4132, I4132
336      symrel(:,:,2) = genrot(:,:)
337      symrel(:,:,3) = genmmp(:,:)
338      tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
339      symrel(:,:,4) = genmpm(:,:)
340      tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
341      symrel(:,:,5) = genswm(:,:)
342      tnons(:,5)=(/0.75d0,0.25d0,0.25d0/)
343      if(shubnikov==3)symafm(5)=-1
344      nogen=5
345    case (215,216,217)        !PB43m, FB43m, IB43m
346      symrel(:,:,2) = genrot(:,:)
347      symrel(:,:,3) = genmmp(:,:)
348      symrel(:,:,4) = genmpm(:,:)
349      symrel(:,:,5) = genswp(:,:)
350      if(shubnikov==3)symafm(5)=-1
351      nogen=5
352    case (218,219)                !PB43n, FB343c
353      symrel(:,:,2) = genrot(:,:)
354      symrel(:,:,3) = genmmp(:,:)
355      symrel(:,:,4) = genmpm(:,:)
356      symrel(:,:,5) = genswp(:,:)
357      tnons(:,5)=(/0.5d0,0.5d0,0.5d0/)
358      if(shubnikov==3)symafm(5)=-1
359      nogen=5
360    case (220)                !IB43d
361      symrel(:,:,2) = genrot(:,:)
362      symrel(:,:,3) = genmmp(:,:)
363      tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
364      symrel(:,:,4) = genmpm(:,:)
365      tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
366      symrel(:,:,5) = genswp(:,:)
367      tnons(:,5)=(/0.25d0,0.25d0,0.25d0/)
368      if(shubnikov==3)symafm(5)=-1
369      nogen=5
370    case (221,225,229)        !Pm3m,Fm3m,Im3m
371      symrel(:,:,2) = genrot(:,:)
372      symrel(:,:,3) = genmmp(:,:)
373      symrel(:,:,4) = genmpm(:,:)
374      symrel(:,:,5) = genswm(:,:)
375      if(shubnikov==3)then
376        if(spgroupma==94 .or. spgroupma==95 .or. spgroupma==118 .or. &
377 &       spgroupma==119 .or. spgroupma==142 .or. spgroupma==143   )symafm(5)=-1
378      end if
379      nogen=5
380    case (222)                !PnB3n
381      if (spgorig==1) then
382        symrel(:,:,2) = genrot(:,:)
383        symrel(:,:,3) = genmmp(:,:)
384        tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
385        symrel(:,:,4) = genmpm(:,:)
386        tnons(:,4)=(/0.5d0,0.d0,0.5d0/)
387        symrel(:,:,5) = genswm(:,:)
388        tnons(:,5)=(/0.d0,0.d0,0.5d0/)
389        if(shubnikov==3)then
390          if(spgroupma==100 .or. spgroupma==101)symafm(5)=-1
391        end if
392        nogen=5
393        nsym=24
394        call bldgrp(msym,nogen,nsym,symafm,symrel,tnons)
395        do ii=1,24
396          symrel(:,:,ii+24)=-symrel(:,:,ii)
397          tnons(:,ii+24)=tnons(:,ii)
398          if(shubnikov==3)then
399            if(spgroupma==100 .or. spgroupma==102)symafm(ii+24)=-symafm(ii)
400            if(spgroupma==101)symafm(ii+24)=symafm(ii)
401          end if
402        end do
403        nogen=0 ; nsym=48
404      else
405        symrel(:,:,2) = genrot(:,:)
406        symrel(:,:,3) = genmmp(:,:)
407        symrel(:,:,4) = genmpm(:,:)
408        symrel(:,:,5) = genswm(:,:)
409        if(shubnikov==3)then
410          if(spgroupma==100 .or. spgroupma==101)symafm(5)=-1
411        end if
412        nogen=5
413        nsym=24
414        call bldgrp(msym,nogen,nsym,symafm,symrel,tnons)
415        do ii=1,24
416          symrel(:,:,ii+24)=-symrel(:,:,ii)
417          tnons(:,ii+24)=(/0.5d0,0.5d0,0.5d0/)
418          if(shubnikov==3)then
419            if(spgroupma==100 .or. spgroupma==102)symafm(ii+24)=-symafm(ii)
420            if(spgroupma==101)symafm(ii+24)=symafm(ii)
421          end if
422        end do
423        nogen=0 ; nsym=48
424      end if
425    case (223,226)           ! PmB3n, FmB3c
426      symrel(:,:,2) = genrot(:,:)
427      symrel(:,:,3) = genmmp(:,:)
428      symrel(:,:,4) = genmpm(:,:)
429      symrel(:,:,5) = genswm(:,:)
430      tnons(:,5)=(/0.5d0,0.5d0,0.5d0/)
431      if(shubnikov==3)then
432        if(spgroupma==106 .or. spgroupma==124 .or. &
433 &       spgroupma==107 .or. spgroupma==125     )symafm(5)=-1
434      end if
435      nogen=5
436    case (224)                !PnB3m
437      if (spgorig==1) then
438        symrel(:,:,2) = genrot(:,:)
439        symrel(:,:,3) = genmmp(:,:)
440        symrel(:,:,4) = genmpm(:,:)
441        symrel(:,:,5) = genswm(:,:)
442        tnons(:,5)=(/0.5d0,0.5d0,0.5d0/)
443        symrel(:,:,6) = genmmm(:,:)
444        tnons(:,6)=(/0.5d0,0.5d0,0.5d0/)
445      else
446        symrel(:,:,2) = genrot(:,:)
447        symrel(:,:,3) = genmmp(:,:)
448        tnons(:,3)=(/0.5d0,0.5d0,0.0d0/)
449        symrel(:,:,4) = genmpm(:,:)
450        tnons(:,4)=(/0.5d0,0.0d0,0.5d0/)
451        symrel(:,:,5) = genswm(:,:)
452        tnons(:,5)=(/0.5d0,0.5d0,0.d0/)
453        symrel(:,:,6) = genmmm(:,:)
454      end if
455      if(shubnikov==3)then
456        if(spgroupma==112 .or. spgroupma==113)symafm(5)=-1
457        if(spgroupma==112 .or. spgroupma==114)symafm(6)=-1
458      end if
459      nogen=6
460    case (227)                !FdB3m
461      if (spgorig==1) then
462        symrel(:,:,2) = genrot(:,:)
463        symrel(:,:,3) = genmmp(:,:)
464        tnons(:,3)=(/0.d0,0.5d0,0.5d0/)
465        symrel(:,:,4) = genmpm(:,:)
466        tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
467        symrel(:,:,5) = genswm(:,:)
468        tnons(:,5)=(/0.75d0,0.25d0,0.75d0/)
469        symrel(:,:,6) = genmmm(:,:)
470        tnons(:,6)=(/0.25d0,0.25d0,0.25d0/)
471      else
472        symrel(:,:,2) = genrot(:,:)
473        symrel(:,:,3) = genmmp(:,:)
474        tnons(:,3)=(/0.75d0,0.25d0,0.5d0/)
475        symrel(:,:,4) = genmpm(:,:)
476        tnons(:,4)=(/0.25d0,0.5d0,0.75d0/)
477        symrel(:,:,5) = genswm(:,:)
478        tnons(:,5)=(/0.75d0,0.25d0,0.5d0/)
479        symrel(:,:,6) = genmmm(:,:)
480      end if
481      if(shubnikov==3)then
482        if(spgroupma==130 .or. spgroupma==131)symafm(5)=-1
483        if(spgroupma==130 .or. spgroupma==132)symafm(6)=-1
484      end if
485      nogen=6
486    case (228)                !FdB3c
487      if (spgorig==1) then
488        symrel(:,:,2) = genrot(:,:)
489        symrel(:,:,3) = genmmp(:,:)
490        tnons(:,3)=(/0.d0,0.5d0,0.5d0/)
491        symrel(:,:,4) = genmpm(:,:)
492        tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
493        symrel(:,:,5) = genswm(:,:)
494        tnons(:,5)=(/0.75d0,0.25d0,0.75d0/)
495        symrel(:,:,6) = genmmm(:,:)
496        tnons(:,6)=(/0.75d0,0.75d0,0.75d0/)
497      else
498        symrel(:,:,2) = genrot(:,:)
499        symrel(:,:,3) = genmmp(:,:)
500        tnons(:,3)=(/0.25d0,0.75d0,0.5d0/)
501        symrel(:,:,4) = genmpm(:,:)
502        tnons(:,4)=(/0.75d0,0.5d0,0.25d0/)
503        symrel(:,:,5) = genswm(:,:)
504        tnons(:,5)=(/0.75d0,0.25d0,0.d0/)
505        symrel(:,:,6) = genmmm(:,:)
506      end if
507      if(shubnikov==3)then
508        if(spgroupma==136 .or. spgroupma==137)symafm(5)=-1
509        if(spgroupma==136 .or. spgroupma==138)symafm(6)=-1
510      end if
511      nogen=6
512    case (230)                !IaB3d
513      symrel(:,:,2) = genrot(:,:)
514      symrel(:,:,3) = genmmp(:,:)
515      tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
516      symrel(:,:,4) = genmpm(:,:)
517      tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
518      symrel(:,:,5) = genswm(:,:)
519      tnons(:,5)=(/0.75d0,0.25d0,0.25d0/)
520      if(shubnikov==3)then
521        if(spgroupma==147 .or. spgroupma==148)symafm(5)=-1
522      end if
523      nogen=5
524    end select
525 
526 !  End CUBIC space groups
527  end if
528 
529 !***************************************************************************
530 
531  call timab(47,1,tsec)
532 
533  if (nogen>0) then
534    call bldgrp(msym,nogen,nsym,symafm,symrel,tnons)
535  end if
536 
537  call timab(47,2,tsec)
538 
539  call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig)
540 
541 end subroutine symsgcube

m_symsg/symsghexa [ Functions ]

[ Top ] [ m_symsg ] [ Functions ]

NAME

 symsghexa

FUNCTION

 Yields all the TRIGONAL & HEXAGONAL symmetry operations starting from the space group symbol.
 according to the International Tables of Crystallography, 1983.

INPUTS

 msym = default number of symmetries
 nsym = the number of symmetry operations
 shubnikov= magnetic type of the space group to be generated
 spgaxor = ossible orientation of the axes system
 spgorig = the origin choice (1 or 2) for the axes system
 spgroup = the numeric symbol of the space groups
 spgroupma= number of the magnetic space group

OUTPUT

 brvltt = bravais lattice type, here, only for rhombohedral groups
  with hexagonal axes
 symafm(nsym)=(anti)ferromagnetic part of symmetry operations
 symrel(3,3,nsym) = 3D matrix containg symmetry operations
 tnons(3,nsym) = 2D matrix containing translations associated

PARENTS

      gensymspgr

CHILDREN

      bldgrp,spgdata

SOURCE

576 subroutine symsghexa(brvltt,msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma,symafm,symrel,tnons)
577 
578 
579 !This section has been created automatically by the script Abilint (TD).
580 !Do not modify the following lines by hand.
581 #undef ABI_FUNC
582 #define ABI_FUNC 'symsghexa'
583 !End of the abilint section
584 
585  implicit none
586 
587 !Arguments ------------------------------------
588 !scalars
589  integer,intent(in) :: msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma
590  integer,intent(inout) :: brvltt !vz_i
591 !arrays
592  integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i
593  real(dp),intent(inout) :: tnons(3,msym) !vz_i
594 
595 !Local variables -----------------------------
596 !scalars
597  integer :: ii,nogen,sporder
598  real(dp),parameter :: fivesixth=5.0d0/6.0d0,twothird=2.0d0/3.0d0
599  character(len=1) :: brvsb
600  character(len=15) :: intsb,ptintsb,ptschsb,schsb
601  character(len=35) :: intsbl
602 !arrays
603  integer :: genm(3,3),genmmp(3,3),genswm(3,3),genswmmm(3,3),genswmmp(3,3)
604  integer :: genswp(3,3)
605 
606 !*************************************************************************
607 
608 !DEBUG
609 !write(std_out,*) 'symsghexa',spgroup,shubnikov,spgroupma
610 !ENDDEBUG
611 
612 !The identity operation belongs to all space groups
613  symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1
614 
615 !Predefine some generators
616  genswm(:,:)=0 ; genswm(2,1)=1 ; genswm(1,2)=1 ; genswm(3,3)=-1 !reshape((/0,1,0,1,0,0,0,0,-1/),(/3,3/),(/0,0/),(/2,1/) )
617  genswmmm(:,:)=0 ; genswmmm(2,1)=-1 ; genswmmm(1,2)=-1 ; genswmmm(3,3)=-1
618 !reshape((/0,-1,0,-1,0,0,0,0,-1/),(/3,3/),(/0,0/),(/2,1/) )
619  genswmmp(:,:)=0 ; genswmmp(2,1)=-1 ; genswmmp(1,2)=-1 ; genswmmp(3,3)=1
620 !reshape((/0,-1,0,-1,0,0,0,0,1/),(/3,3/),(/0,0/),(/2,1/) )
621  genswp(:,:)=0 ; genswp(2,1)=1 ; genswp(1,2)=1 ; genswp(3,3)=1
622 !reshape((/0,1,0,1,0,0,0,0,1/),(/3,3/),(/0,0/),(/2,1/) )
623  genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)=1
624 !reshape((/-1,0,0,0,-1,0,0,0,1/),(/3,3/),(/0,0/),(/2,1/) )
625 
626 !Initialize the associated translations matrix to 0
627  do ii=1,nsym
628    tnons(:,ii)= 0.0d0
629  end do
630  nogen=0
631 
632 !Default non-magnetic behaviour
633  symafm(1:nsym)=1
634 
635 !*************************************************************************
636 
637 !Treat TRIGONAL case
638  if(143<=spgroup .and. spgroup<=167)then
639 
640 !  The hexagonal axis choice (orientation) is first treated
641    if (spgaxor == 1) then
642 
643 !    This matrix is common to ALL trigonal spatial groups in this orientation
644 !    (Note : this is the 3- symmetry operation)
645      symrel(:,:,2)=0 ; symrel(1,1,2)=-1 ; symrel(1,2,2)=1 ; symrel(2,1,2)=-1 ; symrel(3,3,2)=1
646 !    reshape((/-1,1,0,-1,0,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) )
647 
648 !    Assigns the generators to each space group
649      select case (spgroup)
650      case (143,146,147,148)        !P3, R3, PB3, RB3
651        symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(2,2,3)=-1 ; symrel(3,3,3)=1
652 !        reshape((/0,-1,0,1,-1,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) )
653        nogen=0
654      case (144)                !P31
655        tnons(:,2)=(/0.d0,0.d0,twothird/)
656        symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(2,2,3)=-1 ; symrel(3,3,3)=1
657 !        reshape((/0,-1,0,1,-1,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) )
658        tnons(:,3)=(/0.d0,0.d0,third/)
659        nogen=0
660      case (145)                !P32
661        tnons(:,2)=(/0.d0,0.d0,third/)
662        symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(2,2,3)=-1 ; symrel(3,3,3)=1
663 !        reshape((/0,-1,0,1,-1,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) )
664        tnons(:,3)=(/0.d0,0.d0,twothird/)
665        nogen=0
666      case (149)                !P312
667        symrel(:,:,3) = genswmmm(:,:)
668        nogen=3
669      case (150,155)                !P321, R32
670        symrel(:,:,3) = genswm(:,:)
671        nogen=3
672      case (151)                !P3112
673        tnons(:,2)=(/0.d0,0.d0,twothird/)
674        symrel(:,:,3) = genswmmm(:,:)
675        nogen=3
676      case (152)                !P3121
677        tnons(:,2)=(/0.d0,0.d0,twothird/)
678        symrel(:,:,3) = genswm(:,:)
679        nogen=3
680      case (153)                !P3212
681        tnons(:,2)=(/0.d0,0.d0,third/)
682        symrel(:,:,3) = genswmmm(:,:)
683        nogen=3
684      case (154)                !P3221
685        tnons(:,2)=(/0.d0,0.d0,third/)
686        symrel(:,:,3) = genswm(:,:)
687        nogen=3
688      case (156,160,164,166)        !P3m1, R3m, PB3m1, RB3m
689        symrel(:,:,3) = genswmmp(:,:)
690        nogen=3
691      case (157,162)                !P31m, PB31m
692        symrel(:,:,3) = genswp(:,:)
693        nogen=3
694      case (158,161,165,167)        !P3c1, R3c, PB3c1, RB3c
695        symrel(:,:,3) = genswmmp(:,:)
696        tnons(:,3)=(/0.d0,0.d0,0.5d0/)
697        nogen=3
698      case (159,163)                !P31c, PB31c
699        symrel(:,:,3) = genswp(:,:)
700        tnons(:,3)=(/0.d0,0.d0,0.5d0/)
701        nogen=3
702      end select
703 
704      select case (spgroup)
705      case (146,148,155,160,166,167)
706        brvltt=7
707      end select
708 
709 !    Quite simple, because the generator of even order is always the third one.
710      if(shubnikov==3)then
711        select case(spgroupma)
712        case (23,27,31,35,39,43,47,51,55,59,63,67,71,76,77,82,83,88,89,94,95,&
713 &         100,101,106,107)
714          symafm(3)=-1
715        end select
716      end if
717 
718    else if (spgaxor == 2) then
719 !    The rhombohedral axis choice (orientation) is now treated
720      write(std_out,*)'rhombohedral axes'
721 !    Assignment of common three-fold rotation
722      symrel(:,:,2)=0 ; symrel(1,3,2)=1 ; symrel(3,2,2)=1 ; symrel(2,1,2)=1
723 !    reshape((/0,0,1,1,0,0,0,1,0/),(/3,3/),(/0,0/),(/2,1/) )
724      symrel(:,:,3)=0 ; symrel(3,1,3)=1 ; symrel(2,3,3)=1 ; symrel(1,2,3)=1
725 !    reshape((/0,1,0,0,0,1,1,0,0/), (/3,3/), (/0,0/), (/2,1/) )
726 
727      select case (spgroup)
728      case (146,148)       !R3
729      case (155,166)       !R32, RB3m
730        symrel(:,:,4) = genswmmm(:,:)
731        nogen=4
732      case (160)           !R3m
733        symrel(:,:,4) = genswp(:,:)
734        nogen=4
735      case (161,167)       !R3c, RB3c
736        symrel(:,:,4) = genswp(:,:)
737        tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
738        nogen=4
739      end select
740 
741      if(shubnikov==3)then
742        select case(spgroupma)
743        case (47,67,71,99,101,106,107)
744          symafm(4)=-1
745        end select
746      end if
747 
748 !    End selection of axis orientation
749    end if
750 
751 !  End trigonal groups
752  end if
753 
754 !*************************************************************************
755 
756 !Treat HEXAGONAL case
757  if(168<=spgroup .and. spgroup<=194)then
758 
759 !  This matrix (6) is common to most hexagonal spatial groups, except 174,187,188,189,190
760    symrel(:,:,2)=0 ; symrel(1,1,2)=1 ; symrel(1,2,2)=-1 ; symrel(2,1,2)=1 ; symrel(3,3,2)=1
761 !  reshape((/1,-1,0,1,0,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) )
762 !  This one (6 bar) is present in the other cases
763    genm(:,:)=0 ; genm(1,2)=-1 ; genm(2,1)=1 ; genm(2,2)=-1 ; genm(3,3)=-1
764 !  reshape((/0,-1,0,1,-1,0,0,0,-1/), (/3,3/), (/0,0/), (/2,1/) )
765    select case(spgroup)
766    case (168,175)        !P6
767      nogen=2
768    case (169)                !P61
769      tnons(:,2)=(/0.d0,0.d0,sixth/)
770      nogen=2
771    case (170)                !P65
772      tnons(:,2)=(/0.d0,0.d0,fivesixth/)
773      nogen=2
774    case (171)                !P62
775      tnons(:,2)=(/0.d0,0.d0,third/)
776      nogen=2
777    case (172)                !P64
778      tnons(:,2)=(/0.d0,0.d0,twothird/)
779      nogen=2
780    case (173,176)                !P63, P63/m
781      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
782      nogen=2
783    case (174)                !PB6
784      symrel(:,:,2) = genm(:,:)
785      nogen=2
786    case (177)                !P622
787      symrel(:,:,3) =  genswm(:,:)
788      nogen=3
789    case (178)                !P6122
790      tnons(:,2)=(/0.d0,0.d0,sixth/)
791      symrel(:,:,3) =  genswm(:,:)
792      nogen=3
793    case (179)                !P6522
794      tnons(:,2)=(/0.d0,0.d0,fivesixth/)
795      symrel(:,:,3) =  genswm(:,:)
796      nogen=3
797    case (180)                !P6222
798      tnons(:,2)=(/0.d0,0.d0,third/)
799      symrel(:,:,3) =  genswm(:,:)
800      nogen=3
801    case (181)                !P6422
802      tnons(:,2)=(/0.d0,0.d0,twothird/)
803      symrel(:,:,3) =  genswm(:,:)
804      nogen=3
805    case (182)                !P6322
806      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
807      symrel(:,:,3) =  genswm(:,:)
808      nogen=3
809    case (183,191)                !P6mm, P6/mmm
810      symrel(:,:,3) = genswp(:,:)
811      nogen=3
812    case (184,192)                !P6cc, P6/mcc
813      symrel(:,:,3) = genswp(:,:)
814      tnons(:,3)=(/0.d0,0.d0,0.5d0/)
815      nogen=3
816    case (185,193)                !P63cm, P63/mcm
817      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
818      symrel(:,:,3) = genswp(:,:)
819      nogen=3
820    case (186,194)                !P63mc, P63/mmc
821      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
822      symrel(:,:,3) = genswp(:,:)
823      tnons(:,3)=(/0.d0,0.d0,0.5d0/)
824      nogen=3
825    case (187)                !PB6m2
826      symrel(:,:,2)=0 ; symrel(1,2,2)=-1 ; symrel(2,1,2)=1 ; symrel(2,2,2)=-1 ; symrel(3,3,2)=-1
827 !      reshape((/0,-1,0,1,-1,0,0,0,-1/), (/3,3/), (/0,0/), (/2,1/) )
828      symrel(:,:,3)=0 ; symrel(1,1,3)=-1 ; symrel(1,2,3)=1 ; symrel(2,2,3)=1 ; symrel(3,3,3)=1
829 !      reshape((/-1,1,0,0,1,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) )
830      nogen=3
831      if (shubnikov==3) then
832        if (spgroupma==211) symafm(2:3)=-1
833        if (spgroupma==212) symafm(2)=-1
834        if (spgroupma==213) symafm(3)=-1
835      end if
836    case (188)                !PB6c2
837      symrel(:,:,2) = genm(:,:)
838      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
839      symrel(:,:,3) = genswmmm(:,:)
840      nogen=3
841    case (189)                !PB62m
842      symrel(:,:,2) = genm(:,:)
843      symrel(:,:,3) = genswp(:,:)
844      nogen=3
845    case (190)                !PB62c
846      symrel(:,:,2) = genm(:,:)
847      symrel(:,:,3) = genswp(:,:)
848      tnons(:,3)=(/0.d0,0.d0,0.5d0/)
849      nogen=3
850    end select
851 
852    if(shubnikov==3)then
853      select case(spgroupma)
854 !      spgroup from 168 to 176 are OK, 177 to 194 are not done
855      case (111,115,119,123,127,131,135,139,141,145,147,152,158,164,170,&
856 &       176,182,187,193,199,205,217,224,230,237,239,247,249,257,259,267,269)
857        symafm(2)=-1
858      case(153,159,165,171,177,183,189,195,&
859 &       201,207,219,225,231,240,241,250,251,260,261,270,271)
860        symafm(3)=-1
861      case(151,157,163,169,175,181,188,194,200,206,218,223,229,236,238,246,248,256,258,266,268)
862        symafm(2:3)=-1
863      end select
864    end if
865 
866    call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig)
867 
868 !  End HEXAGONAL groups
869  end if
870 
871 !***************************************************************************
872 
873 !DEBUG
874 !write(std_out,*) 'symsghexa : out with nogen = ',nogen
875 !ENDDEBUG
876 
877 
878  if (nogen>0) then
879    call bldgrp(msym,nogen,nsym,symafm,symrel,tnons)
880  end if
881 
882 !DEBUG
883 !write(std_out,*)'symrel:'
884 !write(std_out,*) symrel(:,:,1:nsym)
885 !ENDDEBUG
886 
887 end subroutine symsghexa

m_symsg/symsgmono [ Functions ]

[ Top ] [ m_symsg ] [ Functions ]

NAME

 symsgmono

FUNCTION

 Yields all the MONOCLINIC symmetry operations starting from the space group symbol.
 according to the International Tables of Crystallography, 1983.
 It solves also the problem of the axes orientation
 according to the spgaxor

INPUTS

 msym = default number of symmetries
 nsym = the number of symmetry operations
 shubnikov= magnetic type of the space group to be generated
 spgaxor = the orientation choice of the unit cell
 spgorig = possible origin of the axes system
 spgroup = the numeric symbol of the space groups
 spgroupma= number of the magnetic space group

OUTPUT

 brvltt = bravais lattice type, here, only for rhombohedral groups
  with hexagonal axes (1=P; 2=I; 3=F; 4=C; 5=A; 6=B; 7=R)
 symafm(nsym)=(anti)ferromagnetic part of symmetry operations
 symrel(3,3,nsym) = 3D matrix containg symmetry operations
 tnons(3,nsym) = 2D matrix containing translations associated

PARENTS

      gensymspgr

CHILDREN

      spgdata

SOURCE

 924 subroutine symsgmono(brvltt,msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma,symafm,symrel,tnons)
 925 
 926 
 927 !This section has been created automatically by the script Abilint (TD).
 928 !Do not modify the following lines by hand.
 929 #undef ABI_FUNC
 930 #define ABI_FUNC 'symsgmono'
 931 !End of the abilint section
 932 
 933  implicit none
 934 
 935 !Arguments ------------------------------------
 936 !scalars
 937  integer,intent(in) :: msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma
 938  integer,intent(inout) :: brvltt !vz_i
 939 !arrays
 940  integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i
 941  real(dp),intent(inout) :: tnons(3,msym) !vz_i
 942 
 943 !Local variables -----------------------------
 944 !scalars
 945  integer :: sporder
 946  character(len=1) :: brvsb
 947  character(len=15) :: intsb,ptintsb,ptschsb,schsb
 948  character(len=35) :: intsbl
 949 !arrays
 950  integer :: genmmm(3,3),genmmp(3,3),genmpm(3,3),genmpp(3,3),genpmm(3,3)
 951  integer :: genpmp(3,3),genppm(3,3)
 952 
 953 ! *************************************************************************
 954 !the identity operation belonging to all space groups
 955  symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1
 956 
 957 !Predefine some generators
 958  genmpp(:,:)=0 ; genmpp(1,1)=-1 ; genmpp(2,2)= 1 ; genmpp(3,3)= 1
 959  genpmp(:,:)=0 ; genpmp(1,1)= 1 ; genpmp(2,2)=-1 ; genpmp(3,3)= 1
 960  genppm(:,:)=0 ; genppm(1,1)= 1 ; genppm(2,2)= 1 ; genppm(3,3)=-1
 961  genpmm(:,:)=0 ; genpmm(1,1)= 1 ; genpmm(2,2)=-1 ; genpmm(3,3)=-1
 962  genmpm(:,:)=0 ; genmpm(1,1)=-1 ; genmpm(2,2)= 1 ; genmpm(3,3)=-1
 963  genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)= 1
 964  genmmm(:,:)=0 ; genmmm(1,1)=-1 ; genmmm(2,2)=-1 ; genmmm(3,3)=-1
 965 
 966 !Default non-magnetic behaviour
 967  symafm(1:nsym)=1
 968 
 969 !assigns the generators to each space group
 970  select case (spgroup)
 971  case (3)                 ! P2
 972    select case (spgaxor)
 973    case (1)                ! 3:b, P2_b = P2
 974      symrel(:,:,2) = genmpm(:,:)
 975    case (2)                ! 3:a, P2_a = P2
 976      symrel(:,:,2) = genpmm(:,:)
 977    case (3)                ! 3:c, P2_c = P2
 978      symrel(:,:,2) = genmmp(:,:)
 979    end select
 980    if(shubnikov==3)symafm(2)=-1
 981  case (4)                ! P21
 982    select case (spgaxor)
 983    case (1)                ! 3:b, P21_b = P21
 984      symrel(:,:,2) = genmpm(:,:)
 985      tnons(:,2)=(/0.d0,0.5d0,0.d0/)
 986    case (2)                ! 3:a, P21_a = P21
 987      symrel(:,:,2) = genpmm(:,:)
 988      tnons(:,2)=(/0.5d0,0.d0,0.d0/)
 989    case (3)                ! 3:c, P21_c = P21
 990      symrel(:,:,2) = genmmp(:,:)
 991      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
 992    end select
 993    if(shubnikov==3)symafm(2)=-1
 994  case (5)                ! C2
 995    select case (spgaxor)
 996    case (1)                ! 5:b1, C2  = C2
 997      symrel(:,:,2) = genmpm(:,:)
 998      brvltt=4
 999    case (2)                ! 5:a1, B2_a = C2
1000      symrel(:,:,2) = genpmm(:,:)
1001      brvltt=6
1002    case (3)                ! 5:a2, C2_a = C2
1003      symrel(:,:,2) = genpmm(:,:)
1004      brvltt=4
1005    case (4)                ! 5:a3, I2_a = C2
1006      symrel(:,:,2) = genpmm(:,:)
1007      brvltt=2
1008    case (5)                ! 5:b2, A2_b = C2
1009      symrel(:,:,2) = genmpm(:,:)
1010      brvltt=5
1011    case (6)                ! 5:b3, I2_b = C2
1012      symrel(:,:,2) = genmpm(:,:)
1013      brvltt=2
1014    case (7)                ! 5:c1, A2_c = C2
1015      symrel(:,:,2) = genmmp(:,:)
1016      brvltt=5
1017    case (8)                ! 5:c2, B2_c = C2
1018      symrel(:,:,2) = genmmp(:,:)
1019      brvltt=6
1020    case (9)                ! 5:c3, I2_c = C2
1021      symrel(:,:,2) = genmmp(:,:)
1022      brvltt=2
1023    end select
1024    if(shubnikov==3)symafm(2)=-1
1025  case (6)                ! Pm
1026    select case (spgaxor)
1027    case (1)                ! Pm_b = Pm
1028      symrel(:,:,2) = genpmp(:,:)
1029    case (2)                ! Pm_a = Pm
1030      symrel(:,:,2) = genmpp(:,:)
1031    case (3)                ! Pm_c = Pm
1032      symrel(:,:,2) = genppm(:,:)
1033    end select
1034    if(shubnikov==3)symafm(2)=-1
1035  case (7)                ! Pc
1036    select case (spgaxor)
1037    case (1)                ! 7:b1, Pc_b = Pc
1038      symrel(:,:,2) = genpmp(:,:)
1039      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1040    case (2)                ! 7:a1, Pb_a = Pc
1041      symrel(:,:,2) = genmpp(:,:)
1042      tnons(:,2)=(/0.d0,0.5d0,0.d0/)
1043    case (3)                ! 7:a2, Pn_a = Pc
1044      symrel(:,:,2) = genmpp(:,:)
1045      tnons(:,2)=(/0.d0,0.5d0,0.5d0/)
1046    case (4)                ! 7:a3, Pc_a = Pc
1047      symrel(:,:,2) = genmpp(:,:)
1048      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1049    case (5)                ! 7:b2, Pn_b = Pc
1050      symrel(:,:,2) = genpmp(:,:)
1051      tnons(:,2)=(/0.5d0,0.d0,0.5d0/)
1052    case (6)                ! 7:b3, Pa_b = Pc
1053      symrel(:,:,2) = genpmp(:,:)
1054      tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1055    case (7)                ! 7:c1, Pa_c = Pc
1056      symrel(:,:,2) = genppm(:,:)
1057      tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1058    case (8)                ! 7:c2, Pn_c = Pc
1059      symrel(:,:,2) = genppm(:,:)
1060      tnons(:,2)=(/0.5d0,0.5d0,0.d0/)
1061    case (9)                ! 7:c3, Pb_c = Pb
1062      symrel(:,:,2) = genppm(:,:)
1063      tnons(:,2)=(/0.d0,0.5d0,0.d0/)
1064    end select
1065    if(shubnikov==3)symafm(2)=-1
1066  case (8)                ! Cm
1067    select case (spgaxor)
1068    case (1)                ! 8:b1, Cm = Cm
1069      symrel(:,:,2) = genpmp(:,:)
1070      brvltt=4
1071    case (2)                ! 8:a1, Bm_a = Cm
1072      symrel(:,:,2) = genmpp(:,:)
1073      brvltt=6
1074    case (3)                ! 8:a2, Cm_a = Cm
1075      symrel(:,:,2) = genmpp(:,:)
1076      brvltt=4
1077    case (4)                ! 8:a3, Im_a = Cm
1078      symrel(:,:,2) = genmpp(:,:)
1079      brvltt=2
1080    case (5)                ! 8:b2, Am_b = Cm
1081      symrel(:,:,2) = genpmp(:,:)
1082      brvltt=5
1083    case (6)                ! 8:b3, Im_b = Cm
1084      symrel(:,:,2) = genpmp(:,:)
1085      brvltt=2
1086    case (7)                ! 8:c1, Am_c = Cm
1087      symrel(:,:,2) = genppm(:,:)
1088      brvltt=5
1089    case (8)                ! 8:c2, Bm_c = Bm
1090      symrel(:,:,2) = genppm(:,:)
1091      brvltt=6
1092    case (9)                ! 8:c3, Im_c = Cm
1093      symrel(:,:,2) = genppm(:,:)
1094      brvltt=2
1095    end select
1096    if(shubnikov==3)symafm(2)=-1
1097  case (9)                ! Cc
1098    select case (spgaxor)
1099    case (1)                ! 9:b1, Cc_b = Cc
1100      symrel(:,:,2) = genpmp(:,:)
1101      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1102      brvltt=4
1103    case (2)                ! 9:a1, Bb_a = Cc
1104      symrel(:,:,2) = genmpp(:,:)
1105      tnons(:,2)=(/0.d0,0.5d0,0.d0/)
1106      brvltt=6
1107    case (3)                ! 9:a2, Cn_a = Cc
1108      symrel(:,:,2) = genmpp(:,:)
1109      tnons(:,2)=(/0.d0,0.5d0,0.5d0/)
1110      brvltt=4
1111    case (4)                ! 9:a3, Ic_a = Cc
1112      symrel(:,:,2) = genmpp(:,:)
1113      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1114      brvltt=2
1115    case (5)                ! 9:b2, An_b = Cc
1116      symrel(:,:,2) = genpmp(:,:)
1117      tnons(:,2)=(/0.5d0,0.d0,0.5d0/)
1118      brvltt=5
1119    case (6)                ! 9:b3, Ia_b = Cc
1120      symrel(:,:,2) = genpmp(:,:)
1121      tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1122      brvltt=2
1123    case (7)                ! 9:c1, Aa_c = Cc
1124      symrel(:,:,2) = genppm(:,:)
1125      tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1126      brvltt=5
1127    case (8)                ! 9:c2, B(b+c)_c = Cc
1128      symrel(:,:,2) = genppm(:,:)
1129      tnons(:,2)=(/0.d0,0.5d0,0.5d0/)
1130      brvltt=6
1131    case (9)                ! 9:c3, Ib_c = Cc
1132      symrel(:,:,2) = genppm(:,:)
1133      tnons(:,2)=(/0.d0,0.5d0,0.d0/)
1134      brvltt=2
1135    end select
1136    if(shubnikov==3)symafm(2)=-1
1137  case (10)                ! P2/m
1138    select case (spgaxor)
1139    case (1)                ! 10:b, P2/m = P2/m
1140      symrel(:,:,2) = genmpm(:,:)
1141      symrel(:,:,3) = genmmm(:,:)
1142      symrel(:,:,4) = genpmp(:,:)
1143    case (2)                ! 10:a, P2/m_a = P2/m
1144      symrel(:,:,2) = genpmm(:,:)
1145      symrel(:,:,3) = genmmm(:,:)
1146      symrel(:,:,4) = genmpp(:,:)
1147    case (3)                ! 10:c, P2/m_c = P2/m
1148      symrel(:,:,2) = genmmp(:,:)
1149      symrel(:,:,3) = genmmm(:,:)
1150      symrel(:,:,4) = genppm(:,:)
1151    end select
1152    if(shubnikov==3)then
1153      symafm(2:4)=-1 ! Default
1154      if(spgroupma==44)symafm(4)=1
1155      if(spgroupma==45)symafm(2)=1
1156      if(spgroupma==46)symafm(3)=1
1157    end if
1158  case (11)                ! P21/m
1159    select case (spgaxor)
1160    case (1)                ! 11:b, P21/m = P21/m
1161      symrel(:,:,2) = genmpm(:,:)
1162      symrel(:,:,3) = genmmm(:,:)
1163      symrel(:,:,4) = genpmp(:,:)
1164      tnons(:,2)=(/0.d0,0.5d0,0.d0/)
1165      tnons(:,4)=(/0.d0,0.5d0,0.d0/)
1166    case (2)                ! 11:a, P21/m_a = P21/m
1167      symrel(:,:,2) = genpmm(:,:)
1168      symrel(:,:,3) = genmmm(:,:)
1169      symrel(:,:,4) = genmpp(:,:)
1170      tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1171      tnons(:,4)=(/0.5d0,0.d0,0.d0/)
1172    case (3)                ! 11:c, P21/m_c = P21/m
1173      symrel(:,:,2) = genmmp(:,:)
1174      symrel(:,:,3) = genmmm(:,:)
1175      symrel(:,:,4) = genppm(:,:)
1176      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1177      tnons(:,4)=(/0.d0,0.d0,0.5d0/)
1178    end select
1179    if(shubnikov==3)then
1180      symafm(2:4)=-1 ! Default
1181      if(spgroupma==52)symafm(4)=1
1182      if(spgroupma==53)symafm(2)=1
1183      if(spgroupma==54)symafm(3)=1
1184    end if
1185  case (12)                ! C2/m
1186    select case (spgaxor)
1187    case (1)                ! 12:b1, C2/m = C2/m
1188      symrel(:,:,2) = genmpm(:,:)
1189      symrel(:,:,3) = genmmm(:,:)
1190      symrel(:,:,4) = genpmp(:,:)
1191      brvltt=4
1192    case (2)                ! 12:a1, B2/m_a = C2/m
1193      symrel(:,:,2) = genpmm(:,:)
1194      symrel(:,:,3) = genmmm(:,:)
1195      symrel(:,:,4) = genmpp(:,:)
1196      brvltt=6
1197    case (3)                ! 12:a2, C2/m_a = C2/m
1198      symrel(:,:,2) = genpmm(:,:)
1199      symrel(:,:,3) = genmmm(:,:)
1200      symrel(:,:,4) = genmpp(:,:)
1201      brvltt=4
1202    case (4)                ! 12:a3, I2/m_a = C2/m
1203      symrel(:,:,2) = genpmm(:,:)
1204      symrel(:,:,3) = genmmm(:,:)
1205      symrel(:,:,4) = genmpp(:,:)
1206      brvltt=2
1207    case (5)                ! 12:b2, A2/m_b = C2/m
1208      symrel(:,:,2) = genmpm(:,:)
1209      symrel(:,:,3) = genmmm(:,:)
1210      symrel(:,:,4) = genpmp(:,:)
1211      brvltt=5
1212    case (6)                ! 12:b3, I2/m_b = C2/m
1213      symrel(:,:,2) = genmpm(:,:)
1214      symrel(:,:,3) = genmmm(:,:)
1215      symrel(:,:,4) = genpmp(:,:)
1216      brvltt=2
1217    case (7)                ! 12:c1, A2/m_c = C2/m
1218      symrel(:,:,2) = genmmp(:,:)
1219      symrel(:,:,3) = genmmm(:,:)
1220      symrel(:,:,4) = genppm(:,:)
1221      brvltt=5
1222    case (8)                ! 12:c2, B2/m_c = B2/m
1223      symrel(:,:,2) = genmmp(:,:)
1224      symrel(:,:,3) = genmmm(:,:)
1225      symrel(:,:,4) = genppm(:,:)
1226      brvltt=6
1227    case (9)                ! 12:c3, I2/m_c = C2/m
1228      symrel(:,:,2) = genmmp(:,:)
1229      symrel(:,:,3) = genmmm(:,:)
1230      symrel(:,:,4) = genppm(:,:)
1231      brvltt=2
1232    end select
1233    if(shubnikov==3)then
1234      symafm(2:4)=-1 ! Default
1235      if(spgroupma==60)symafm(4)=1
1236      if(spgroupma==61)symafm(2)=1
1237      if(spgroupma==62)symafm(3)=1
1238    end if
1239  case (13)                ! P2/c
1240    select case (spgaxor)
1241    case (1)                ! 13:b1, P2/c = P2/c
1242      symrel(:,:,2) = genmpm(:,:)
1243      symrel(:,:,3) = genmmm(:,:)
1244      symrel(:,:,4) = genpmp(:,:)
1245      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1246      tnons(:,4)=(/0.d0,0.d0,0.5d0/)
1247    case (2)                ! 13:a1, P2/b_a = P2/c
1248      symrel(:,:,2) = genpmm(:,:)
1249      symrel(:,:,3) = genmmm(:,:)
1250      symrel(:,:,4) = genmpp(:,:)
1251      tnons(:,2)=(/0.d0,0.5d0,0.d0/)
1252      tnons(:,4)=(/0.d0,0.5d0,0.d0/)
1253    case (3)                ! 13:a2, P2/n_a = P2/c
1254      symrel(:,:,2) = genpmm(:,:)
1255      symrel(:,:,3) = genmmm(:,:)
1256      symrel(:,:,4) = genmpp(:,:)
1257      tnons(:,2)=(/0.d0,0.5d0,0.5d0/)
1258      tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
1259    case (4)                ! 13:a3, P2/c_a = P2/c
1260      symrel(:,:,2) = genpmm(:,:)
1261      symrel(:,:,3) = genmmm(:,:)
1262      symrel(:,:,4) = genmpp(:,:)
1263      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1264      tnons(:,4)=(/0.d0,0.d0,0.5d0/)
1265    case (5)                ! 13:b2, P2/n_b = P2/c
1266      symrel(:,:,2) = genmpm(:,:)
1267      symrel(:,:,3) = genmmm(:,:)
1268      symrel(:,:,4) = genpmp(:,:)
1269      tnons(:,2)=(/0.5d0,0.d0,0.5d0/)
1270      tnons(:,4)=(/0.5d0,0.d0,0.5d0/)
1271    case (6)                ! 13:b3, P2/a_b = P2/c
1272      symrel(:,:,2) = genmpm(:,:)
1273      symrel(:,:,3) = genmmm(:,:)
1274      symrel(:,:,4) = genpmp(:,:)
1275      tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1276      tnons(:,4)=(/0.5d0,0.d0,0.d0/)
1277    case (7)                ! 13:c1, P2/a_c = P2/c
1278      symrel(:,:,2) = genmmp(:,:)
1279      symrel(:,:,3) = genmmm(:,:)
1280      symrel(:,:,4) = genppm(:,:)
1281      tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1282      tnons(:,4)=(/0.5d0,0.d0,0.d0/)
1283    case (8)                ! 13:c2, P2/n_c = P2/c
1284      symrel(:,:,2) = genmmp(:,:)
1285      symrel(:,:,3) = genmmm(:,:)
1286      symrel(:,:,4) = genppm(:,:)
1287      tnons(:,2)=(/0.5d0,0.5d0,0.d0/)
1288      tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
1289    case (9)                ! 13:c3, P2/b_c = P2/b
1290      symrel(:,:,2) = genmmp(:,:)
1291      symrel(:,:,3) = genmmm(:,:)
1292      symrel(:,:,4) = genppm(:,:)
1293      tnons(:,2)=(/0.d0,0.5d0,0.d0/)
1294      tnons(:,4)=(/0.d0,0.5d0,0.d0/)
1295    end select
1296    if(shubnikov==3)then
1297      symafm(2:4)=-1 ! Default
1298      if(spgroupma==67)symafm(4)=1
1299      if(spgroupma==68)symafm(2)=1
1300      if(spgroupma==69)symafm(3)=1
1301    end if
1302  case (14)              ! P21/c
1303    select case (spgaxor)
1304    case (1)             ! 14:b1, P21/c_b = P21/c
1305      symrel(:,:,2) = genmpm(:,:)
1306      symrel(:,:,3) = genmmm(:,:)
1307      symrel(:,:,4) = genpmp(:,:)
1308      tnons(:,2)=(/0.d0,0.5d0,0.5d0/)
1309      tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
1310    case (2)             ! 14:a1, P21/a_b = P21/c
1311      symrel(:,:,2) = genpmm(:,:)
1312      symrel(:,:,3) = genmmm(:,:)
1313      symrel(:,:,4) = genmpp(:,:)
1314      tnons(:,2)=(/0.5d0,0.5d0,0.d0/)
1315      tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
1316    case (3)                ! 14:a2, P21/n_a = P21/c
1317      symrel(:,:,2) = genpmm(:,:)
1318      symrel(:,:,3) = genmmm(:,:)
1319      symrel(:,:,4) = genmpp(:,:)
1320      tnons(:,2)=(/0.5d0,0.5d0,0.5d0/)
1321      tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
1322    case (4)                ! 14:a3, P21/c_a = P21/c
1323      symrel(:,:,2) = genpmm(:,:)
1324      symrel(:,:,3) = genmmm(:,:)
1325      symrel(:,:,4) = genmpp(:,:)
1326      tnons(:,2)=(/0.5d0,0.d0,0.5d0/)
1327      tnons(:,4)=(/0.5d0,0.d0,0.5d0/)
1328    case (5)                ! 14:b2, P21/n_b = P21/c
1329      symrel(:,:,2) = genmpm(:,:)
1330      symrel(:,:,3) = genmmm(:,:)
1331      symrel(:,:,4) = genpmp(:,:)
1332      tnons(:,2)=(/0.5d0,0.5d0,0.5d0/)
1333      tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
1334    case (6)                ! 14:b3, P21/a_b = P21/c
1335      symrel(:,:,2) = genmpm(:,:)
1336      symrel(:,:,3) = genmmm(:,:)
1337      symrel(:,:,4) = genpmp(:,:)
1338      tnons(:,2)=(/0.5d0,0.5d0,0.d0/)
1339      tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
1340    case (7)                ! 14:c1, P21/a_c = P21/c
1341      symrel(:,:,2) = genmmp(:,:)
1342      symrel(:,:,3) = genmmm(:,:)
1343      symrel(:,:,4) = genppm(:,:)
1344      tnons(:,2)=(/0.5d0,0.d0,0.5d0/)
1345      tnons(:,4)=(/0.5d0,0.d0,0.5d0/)
1346    case (8)                ! 14:c2, P21/n_c = P21/c
1347      symrel(:,:,2) = genmmp(:,:)
1348      symrel(:,:,3) = genmmm(:,:)
1349      symrel(:,:,4) = genppm(:,:)
1350      tnons(:,2)=(/0.5d0,0.5d0,0.5d0/)
1351      tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
1352    case (9)                ! 14/c3, P21/b_c = P21/b
1353      symrel(:,:,2) = genmmp(:,:)
1354      symrel(:,:,3) = genmmm(:,:)
1355      symrel(:,:,4) = genppm(:,:)
1356      tnons(:,2)=(/0.d0,0.5d0,0.5d0/)
1357      tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
1358    end select
1359    if(shubnikov==3)then
1360      symafm(2:4)=-1 ! Default
1361      if(spgroupma==77)symafm(4)=1
1362      if(spgroupma==78)symafm(2)=1
1363      if(spgroupma==79)symafm(3)=1
1364    end if
1365  case (15)                ! C2/c
1366    select case (spgaxor)
1367    case (1)                ! 15:b1, C2/c_b = C2/c
1368      symrel(:,:,2) = genmpm(:,:)
1369      symrel(:,:,3) = genmmm(:,:)
1370      symrel(:,:,4) = genpmp(:,:)
1371      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1372      tnons(:,4)=(/0.d0,0.d0,0.5d0/)
1373      brvltt = 4
1374    case (2)                ! 15:a1, B2/b_a = C2/c
1375      symrel(:,:,2) = genpmm(:,:)
1376      symrel(:,:,3) = genmmm(:,:)
1377      symrel(:,:,4) = genmpp(:,:)
1378      tnons(:,2)=(/0.d0,0.5d0,0.d0/)
1379      tnons(:,4)=(/0.d0,0.5d0,0.d0/)
1380      brvltt = 6
1381    case (3)                ! 15:a2, C2/n_a = C2/c
1382      symrel(:,:,2) = genpmm(:,:)
1383      symrel(:,:,3) = genmmm(:,:)
1384      symrel(:,:,4) = genmpp(:,:)
1385      tnons(:,2)=(/0.d0,0.5d0,0.5d0/)
1386      tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
1387      brvltt = 4
1388    case (4)                ! 15:a3, I2/c_a = C2/c
1389      symrel(:,:,2) = genpmm(:,:)
1390      symrel(:,:,3) = genmmm(:,:)
1391      symrel(:,:,4) = genmpp(:,:)
1392      tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1393      tnons(:,4)=(/0.d0,0.d0,0.5d0/)
1394      brvltt = 2
1395    case (5)                ! 15:b2, A2/n_b = C2/c
1396      symrel(:,:,2) = genmpm(:,:)
1397      symrel(:,:,3) = genmmm(:,:)
1398      symrel(:,:,4) = genpmp(:,:)
1399      tnons(:,2)=(/0.5d0,0.d0,0.5d0/)
1400      tnons(:,4)=(/0.5d0,0.d0,0.5d0/)
1401      brvltt = 5
1402    case (6)                ! 15:b3, I2/a_b = C2/c
1403      symrel(:,:,2) = genmpm(:,:)
1404      symrel(:,:,3) = genmmm(:,:)
1405      symrel(:,:,4) = genpmp(:,:)
1406      tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1407      tnons(:,4)=(/0.5d0,0.d0,0.d0/)
1408      brvltt = 2
1409    case (7)                ! 15:c1, A2/a_c = C2/c
1410      symrel(:,:,2) = genmmp(:,:)
1411      symrel(:,:,3) = genmmm(:,:)
1412      symrel(:,:,4) = genppm(:,:)
1413      tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1414      tnons(:,4)=(/0.5d0,0.d0,0.d0/)
1415      brvltt = 5
1416    case (8)                ! 15:c2, B21/b_c = C2/c
1417      symrel(:,:,2) = genmmp(:,:)
1418      symrel(:,:,3) = genmmm(:,:)
1419      symrel(:,:,4) = genppm(:,:)
1420      tnons(:,2)=(/0.d0,0.5d0,0.5d0/)
1421      tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
1422      brvltt = 6
1423    case (9)                ! 15:c3, I2/b_c = C2/c
1424      symrel(:,:,2) = genmmp(:,:)
1425      symrel(:,:,3) = genmmm(:,:)
1426      symrel(:,:,4) = genppm(:,:)
1427      tnons(:,2)=(/0.d0,0.5d0,0.d0/)
1428      tnons(:,4)=(/0.d0,0.5d0,0.d0/)
1429      brvltt = 2
1430    end select
1431    if(shubnikov==3)then
1432      symafm(2:4)=-1 ! Default
1433      if(spgroupma==87)symafm(4)=1
1434      if(spgroupma==88)symafm(2)=1
1435      if(spgroupma==89)symafm(3)=1
1436    end if
1437  end select
1438 
1439  call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig)
1440 
1441 end subroutine symsgmono

m_symsg/symsgortho [ Functions ]

[ Top ] [ m_symsg ] [ Functions ]

NAME

 symsgortho

FUNCTION

 Yields all the ORTHORHOMBIC symmetry operations starting from the space group symbol.
 It deals only with the orthorhombic groups
 taken in the standard orientation
 according to the International Tables of Crystallography, 1983.

INPUTS

 msym = default number of symmetries
 nsym = the number of symmetry operations
 shubnikov= magnetic type of the space group to be generated
 spgorig = the origin choice (1 or 2) for the axes system
 spgaxor = the possible orientation of the axes system
 spgroup = the numeric symbol of the space groups
 spgroupma= number of the magnetic space group

OUTPUT

 symafm(nsym)=(anti)ferromagnetic part of symmetry operations
 symrel(3,3,nsym) = 3D matrix containg symmetry operations
 tnons(3,nsym) = 2D matrix containing translations associated

PARENTS

      gensymspgr

CHILDREN

      bldgrp,spgdata

SOURCE

1476 subroutine symsgortho(msym,nsym,shubnikov,spgaxor,spgorig,spgroup,&
1477 &   spgroupma,symafm,symrel,tnons)
1478 
1479 
1480 !This section has been created automatically by the script Abilint (TD).
1481 !Do not modify the following lines by hand.
1482 #undef ABI_FUNC
1483 #define ABI_FUNC 'symsgortho'
1484 !End of the abilint section
1485 
1486  implicit none
1487 
1488 !Arguments ------------------------------------
1489 !scalars
1490  integer,intent(in) :: msym,nsym,shubnikov,spgaxor,spgorig,spgroup
1491  integer,intent(in) :: spgroupma
1492 !arrays
1493  integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i
1494  real(dp),intent(inout) :: tnons(3,msym) !vz_i
1495 
1496 !Local variables ------------------------------
1497 ! nogen = number of generators selected
1498 !scalars
1499  integer :: nogen,sporder
1500  character(len=1) :: brvsb
1501  character(len=15) :: intsb,ptintsb,ptschsb,schsb
1502  character(len=35) :: intsbl
1503 !arrays
1504  integer :: genmmm(3,3),genmmp(3,3),genmpm(3,3),genmpp(3,3),genpmm(3,3)
1505  integer :: genpmp(3,3),genppm(3,3)
1506 
1507 ! *************************************************************************
1508 
1509 !DEBUG
1510 !write(std_out,*)'symsgortho ( orthorhombic groups) : enter with space group ',spgroup
1511 !ENDDEBUG
1512 
1513 !The orientation of the space group:
1514 !first we will permute the input coordinates of the atoms, xred
1515 !then we will make the calculation in the "normal" space group
1516 !then the coordinates are reoriented to match the initial orientation
1517 !and finally the symrel is reoriented to correspond to the new orientation
1518 !further all the calculations are performed into the space group
1519 !with the user-defined orientation
1520 
1521 !The identity operation belongs to all space groups
1522  symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1
1523 
1524  nogen=4
1525 
1526 !Predefine some generators
1527  genmpp(:,:)=0 ; genmpp(1,1)=-1 ; genmpp(2,2)= 1 ; genmpp(3,3)= 1
1528  genpmp(:,:)=0 ; genpmp(1,1)= 1 ; genpmp(2,2)=-1 ; genpmp(3,3)= 1
1529  genppm(:,:)=0 ; genppm(1,1)= 1 ; genppm(2,2)= 1 ; genppm(3,3)=-1
1530  genpmm(:,:)=0 ; genpmm(1,1)= 1 ; genpmm(2,2)=-1 ; genpmm(3,3)=-1
1531  genmpm(:,:)=0 ; genmpm(1,1)=-1 ; genmpm(2,2)= 1 ; genmpm(3,3)=-1
1532  genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)= 1
1533  genmmm(:,:)=0 ; genmmm(1,1)=-1 ; genmmm(2,2)=-1 ; genmmm(3,3)=-1
1534 
1535 
1536 !For all the groups in this routine symrel(:,:,2) is the same
1537  symrel(:,:,2)=genmmp(:,:)
1538 
1539 !Default non-magnetic behaviour
1540  symafm(1:nsym)=1
1541 
1542 !DEBUG
1543 !write(std_out,*) 'symsgortho:',spgroup,shubnikov,spgroupma
1544 !ENDDEBUG
1545 
1546 !assigns the generators to each space group
1547  select case (spgroup)
1548 !  ORTHORHOMBIC space groups
1549  case (16,21,22,23)        !P222, C222, F222, I222
1550    symrel(:,:,3) = genmpm(:,:)
1551    symrel(:,:,4) = genpmm(:,:)
1552    if(shubnikov==3)symafm(3:4)=-1
1553  case (17,20)                !P2221, C2221
1554    tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1555    symrel(:,:,3) = genpmm(:,:)
1556    symrel(:,:,4) = genmpm(:,:)
1557    tnons(:,4)=(/0.d0,0.d0,0.5d0/)
1558    if(shubnikov==3)then
1559      symafm(4)=-1
1560      if(spgroupma==9) symafm(3)=-1
1561      if(spgroupma==10)symafm(2)=-1
1562      if(spgroupma==33)symafm(3:4)=-1
1563      if(spgroupma==34)symafm(2)=-1
1564    end if
1565  case (18)                !P21212
1566    symrel(:,:,3) = genmpm(:,:)
1567    tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
1568    symrel(:,:,4) = genpmm(:,:)
1569    tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
1570    if(shubnikov==3)then
1571      symafm(3)=-1
1572      if(spgroupma==18)symafm(4)=-1
1573      if(spgroupma==19)symafm(2)=-1
1574    end if
1575  case (19,24)                !P212121, I212121
1576    tnons(:,2)=(/0.5d0,0.d0,0.5d0/)
1577    symrel(:,:,3) = genmpm(:,:)
1578    tnons(:,3)=(/0.d0,0.5d0,0.5d0/)
1579    symrel(:,:,4) = genpmm(:,:)
1580    tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
1581    if(shubnikov==3)symafm(3:4)=-1
1582  case (25,35,38,42,44)        !Pmm2, Cmm2, Amm2, Fmm2, Imm2
1583    symrel(:,:,3) = genpmp(:,:)
1584    symrel(:,:,4) = genmpp(:,:)
1585  case (26,36)                !Pmc21, Cmc21
1586    tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1587    symrel(:,:,3) = genpmp(:,:)
1588    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
1589    symrel(:,:,4) = genmpp(:,:)
1590  case (27,37)                !Pcc2, Ccc2
1591    symrel(:,:,3) = genpmp(:,:)
1592    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
1593    symrel(:,:,4) = genmpp(:,:)
1594    tnons(:,4)=(/0.d0,0.d0,0.5d0/)
1595  case (28,40,46)        !Pma2, Ama2, Ima2
1596    symrel(:,:,3) = genpmp(:,:)
1597    tnons(:,3)=(/0.5d0,0.d0,0.d0/)
1598    symrel(:,:,4) = genmpp(:,:)
1599    tnons(:,4)=(/0.5d0,0.d0,0.d0/)
1600  case (29)                 !Pca21
1601    tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1602    symrel(:,:,3) = genpmp(:,:)
1603    tnons(:,3)=(/0.5d0,0.d0,0.d0/)
1604    symrel(:,:,4) = genmpp(:,:)
1605    tnons(:,4)=(/0.5d0,0.d0,0.5d0/)
1606  case (30)                !Pnc2
1607    symrel(:,:,3) = genpmp(:,:)
1608    tnons(:,3)=(/0.d0,0.5d0,0.5d0/)
1609    symrel(:,:,4) = genmpp(:,:)
1610    tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
1611  case (31)                !Pmn21
1612    tnons(:,2)=(/0.5d0,0.d0,0.5d0/)
1613    symrel(:,:,3) = genpmp(:,:)
1614    tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
1615    symrel(:,:,4) = genmpp(:,:)
1616  case (32,41,45)        !Pba2, Aba2, Iba2
1617    symrel(:,:,3) = genpmp(:,:)
1618    tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
1619    symrel(:,:,4) = genmpp(:,:)
1620    tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
1621  case (33)                !Pna21
1622    tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1623    symrel(:,:,3) = genpmp(:,:)
1624    tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
1625    symrel(:,:,4) = genmpp(:,:)
1626    tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
1627  case (34)                !Pnn2
1628    symrel(:,:,3) = genpmp(:,:)
1629    tnons(:,3)=(/0.5d0,0.5d0,0.5d0/)
1630    symrel(:,:,4) = genmpp(:,:)
1631    tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
1632  case (39)                !Abm2
1633    symrel(:,:,3) = genpmp(:,:)
1634    tnons(:,3)=(/0.d0,0.5d0,0.d0/)
1635    symrel(:,:,4) = genmpp(:,:)
1636    tnons(:,4)=(/0.d0,0.5d0,0.d0/)
1637  case (43)                !Fdd2
1638    symrel(:,:,3) = genpmp(:,:)
1639    tnons(:,3)=(/0.25d0,0.25d0,0.25d0/)
1640    symrel(:,:,4) = genmpp(:,:)
1641    tnons(:,4)=(/0.25d0,0.25d0,0.25d0/)
1642  case (47,65,69,71)        !Pmmm, Cmmm, Fmmm, Immm
1643    symrel(:,:,3) = genmpm(:,:)
1644    symrel(:,:,4) = genpmm(:,:)
1645  case (48)                !Pnnn
1646    symrel(:,:,3) = genmpm(:,:)
1647    symrel(:,:,4) = genpmm(:,:)
1648    symrel(:,:,5) = genmmm(:,:)
1649    symrel(:,:,6) = genppm(:,:)
1650    symrel(:,:,7) = genpmp(:,:)
1651    symrel(:,:,8) = genmpp(:,:)
1652    if (spgorig==1) then
1653      tnons(:,5)=(/0.5d0,0.5d0,0.5d0/)
1654      tnons(:,6)=(/0.5d0,0.5d0,0.5d0/)
1655      tnons(:,7)=(/0.5d0,0.5d0,0.5d0/)
1656      tnons(:,8)=(/0.5d0,0.5d0,0.5d0/)
1657    else
1658      tnons(:,2)=(/0.5d0,0.5d0,0.d0/)
1659      tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
1660      tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
1661      tnons(:,6)=(/0.5d0,0.5d0,0.d0/)
1662      tnons(:,7)=(/0.5d0,0.d0,0.5d0/)
1663      tnons(:,8)=(/0.d0,0.5d0,0.5d0/)
1664    end if
1665    if(shubnikov==3)then
1666      if(spgroupma==259)then
1667        symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1
1668      else if(spgroupma==260)then
1669        symafm(3:4)=-1 ; symafm(7:8)=-1
1670      else if(spgroupma==261)then
1671        symafm(5:8)=-1
1672      end if
1673    end if
1674    nogen=0
1675  case (49,66)                !Pccm, Cccm
1676    symrel(:,:,4) = genpmm(:,:)
1677    tnons(:,4)=(/0.d0,0.d0,0.5d0/)
1678    symrel(:,:,3) = genmpm(:,:)
1679    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
1680  case (50)                !Pban
1681    symrel(:,:,3) = genmpm(:,:)
1682    symrel(:,:,4) = genpmm(:,:)
1683    symrel(:,:,5) = genmmm(:,:)
1684    symrel(:,:,6) = genppm(:,:)
1685    symrel(:,:,7) = genpmp(:,:)
1686    symrel(:,:,8) = genmpp(:,:)
1687    if (spgorig==1) then
1688      tnons(:,5)=(/0.5d0,0.5d0,0.d0/)
1689      tnons(:,6)=(/0.5d0,0.5d0,0.d0/)
1690      tnons(:,7)=(/0.5d0,0.5d0,0.d0/)
1691      tnons(:,8)=(/0.5d0,0.5d0,0.d0/)
1692    else
1693      tnons(:,2)=(/0.5d0,0.5d0,0.d0/)
1694      tnons(:,6)=(/0.5d0,0.5d0,0.d0/)
1695      tnons(:,3)=(/0.5d0,0.d0,0.d0/)
1696      tnons(:,7)=(/0.5d0,0.d0,0.d0/)
1697      tnons(:,4)=(/0.d0,0.5d0,0.d0/)
1698      tnons(:,8)=(/0.d0,0.5d0,0.d0/)
1699    end if
1700    if(shubnikov==3)then
1701      if(spgroupma==279)then
1702        symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1
1703      else if(spgroupma==280)then
1704        symafm(3:4)=-1 ; symafm(5:6)=-1
1705      else if(spgroupma==281)then
1706        symafm(3:4)=-1 ; symafm(7:8)=-1
1707      else if(spgroupma==282)then
1708        symafm(2:8:2)=-1
1709      else if(spgroupma==283)then
1710        symafm(5:8)=-1
1711      end if
1712    end if
1713    nogen=0
1714  case (51)                !Pmma
1715    tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1716    symrel(:,:,3) = genmpm(:,:)
1717    symrel(:,:,4) = genpmm(:,:)
1718    tnons(:,4)=(/0.5d0,0.d0,0.d0/)
1719  case (52)                !Pnna
1720    tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1721    symrel(:,:,3) = genmpm(:,:)
1722    tnons(:,3)=(/0.5d0,0.5d0,0.5d0/)
1723    symrel(:,:,4) = genpmm(:,:)
1724    tnons(:,4)=(/0.d0,0.5d0,0.5d0/)
1725  case (53)                 !Pmna
1726    tnons(:,2)=(/0.5d0,0.d0,0.5d0/)
1727    symrel(:,:,3) = genmpm(:,:)
1728    tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
1729    symrel(:,:,4) = genpmm(:,:)
1730  case (54)                !Pcca
1731    tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1732    symrel(:,:,3) = genmpm(:,:)
1733    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
1734    symrel(:,:,4) = genpmm(:,:)
1735    tnons(:,4)=(/0.5d0,0.d0,0.5d0/)
1736  case (55,72)                !Pbam, Ibam
1737    symrel(:,:,3) = genmpm(:,:)
1738    tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
1739    symrel(:,:,4) = genpmm(:,:)
1740    tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
1741  case (56)                !Pccn
1742    tnons(:,2)=(/0.5d0,0.5d0,0.d0/)
1743    symrel(:,:,3) = genmpm(:,:)
1744    tnons(:,3)=(/0.d0,0.5d0,0.5d0/)
1745    symrel(:,:,4) = genpmm(:,:)
1746    tnons(:,4)=(/0.5d0,0.d0,0.5d0/)
1747  case (57)                !Pbcm
1748    tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1749    symrel(:,:,3) = genmpm(:,:)
1750    tnons(:,3)=(/0.d0,0.5d0,0.5d0/)
1751    symrel(:,:,4) = genpmm(:,:)
1752    tnons(:,4)=(/0.d0,0.5d0,0.d0/)
1753  case (58)                !Pnnm
1754    symrel(:,:,3) = genmpm(:,:)
1755    tnons(:,3)=(/0.5d0,0.5d0,0.5d0/)
1756    symrel(:,:,4) = genpmm(:,:)
1757    tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
1758  case (59)                !Pmmn
1759    symrel(:,:,3) = genmpm(:,:)
1760    symrel(:,:,4) = genpmm(:,:)
1761    symrel(:,:,5) = genmmm(:,:)
1762    symrel(:,:,6) = genppm(:,:)
1763    symrel(:,:,7) = genpmp(:,:)
1764    symrel(:,:,8) = genmpp(:,:)
1765    if (spgorig==1) then
1766      tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
1767      tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
1768      tnons(:,5)=(/0.5d0,0.5d0,0.d0/)
1769      tnons(:,6)=(/0.5d0,0.5d0,0.d0/)
1770    else
1771      tnons(:,2)=(/0.5d0,0.5d0,0.d0/)
1772      tnons(:,3)=(/0.d0,0.5d0,0.d0/)
1773      tnons(:,4)=(/0.5d0,0.d0,0.d0/)
1774      tnons(:,6)=(/0.5d0,0.5d0,0.d0/)
1775      tnons(:,7)=(/0.d0,0.5d0,0.d0/)
1776      tnons(:,8)=(/0.5d0,0.d0,0.d0/)
1777    end if
1778    if(shubnikov==3)then
1779      if(spgroupma==407)then
1780        symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1
1781      else if(spgroupma==408)then
1782        symafm(3:4)=-1 ; symafm(5:6)=-1
1783      else if(spgroupma==409)then
1784        symafm(3:4)=-1 ; symafm(7:8)=-1
1785      else if(spgroupma==410)then
1786        symafm(2:8:2)=-1
1787      else if(spgroupma==411)then
1788        symafm(5:8)=-1
1789      end if
1790    end if
1791    nogen=0
1792  case (60)                !Pbcn
1793    tnons(:,2)=(/0.5d0,0.5d0,0.5d0/)
1794    symrel(:,:,3) = genmpm(:,:)
1795    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
1796    symrel(:,:,4) = genpmm(:,:)
1797    tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
1798  case (61,73)                !Pbca, Ibca
1799    tnons(:,2)=(/0.5d0,0.d0,0.5d0/)
1800    symrel(:,:,3) = genmpm(:,:)
1801    tnons(:,3)=(/0.d0,0.5d0,0.5d0/)
1802    symrel(:,:,4) = genpmm(:,:)
1803    tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
1804  case (62)                !Pnma
1805    tnons(:,2)=(/0.5d0,0.d0,0.5d0/)
1806    symrel(:,:,3) = genmpm(:,:)
1807    tnons(:,3)=(/0.d0,0.5d0,0.d0/)
1808    symrel(:,:,4) = genpmm(:,:)
1809    tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
1810  case (63)                !Cmcm
1811    tnons(:,2)=(/0.d0,0.d0,0.5d0/)
1812    symrel(:,:,3) = genmpm(:,:)
1813    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
1814    symrel(:,:,4) = genpmm(:,:)
1815    if(shubnikov==3)then
1816      if(spgroupma==459 .or. spgroupma==463)symafm(2:3)=-1
1817      if(spgroupma==460 .or. spgroupma==464)symafm(4)=-1
1818      if(spgroupma==460 .or. spgroupma==464)symafm(2)=-1
1819      if(spgroupma==461 .or. spgroupma==462)symafm(3:4)=-1
1820    end if
1821  case (64)                !Cmca
1822    tnons(:,2)=(/0.d0,0.5d0,0.5d0/)
1823    symrel(:,:,3) = genmpm(:,:)
1824    tnons(:,3)=(/0.d0,0.5d0,0.5d0/)
1825    symrel(:,:,4) = genpmm(:,:)
1826  case (67)                !Cmma
1827    tnons(:,2)=(/0.d0,0.5d0,0.d0/)
1828    symrel(:,:,4) = genpmm(:,:)
1829    symrel(:,:,3) = genmpm(:,:)
1830    tnons(:,3)=(/0.d0,0.5d0,0.d0/)
1831  case (68)                !Ccca
1832    symrel(:,:,3) = genmpm(:,:)
1833    symrel(:,:,4) = genpmm(:,:)
1834    symrel(:,:,5) = genmmm(:,:)
1835    symrel(:,:,6) = genppm(:,:)
1836    symrel(:,:,7) = genpmp(:,:)
1837    symrel(:,:,8) = genmpp(:,:)
1838    if (spgorig==1) then
1839      tnons(:,2)=(/0.5d0,0.5d0,0.d0/)
1840      tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
1841      tnons(:,5)=(/0.d0,0.5d0,0.5d0/)
1842      tnons(:,6)=(/0.5d0,0.d0,0.5d0/)
1843      tnons(:,7)=(/0.d0,0.5d0,0.5d0/)
1844      tnons(:,8)=(/0.5d0,0.d0,0.5d0/)
1845    else
1846      tnons(:,2)=(/0.5d0,0.d0,0.d0/)
1847      tnons(:,3)=(/0.d0,0.d0,0.5d0/)
1848      tnons(:,4)=(/0.5d0,0.d0,0.5d0/)
1849      tnons(:,6)=(/0.5d0,0.d0,0.d0/)
1850      tnons(:,7)=(/0.d0,0.d0,0.5d0/)
1851      tnons(:,8)=(/0.5d0,0.d0,0.5d0/)
1852    end if
1853    if(shubnikov==3)then
1854      if(spgroupma==513)then
1855        symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1
1856      else if(spgroupma==514)then
1857        symafm(3:4)=-1 ; symafm(5:6)=-1
1858      else if(spgroupma==515)then
1859        symafm(3:4)=-1 ; symafm(7:8)=-1
1860      else if(spgroupma==516)then
1861        symafm(2:8:2)=-1
1862      else if(spgroupma==517)then
1863        symafm(5:8)=-1
1864      end if
1865    end if
1866    nogen=0
1867  case (70)                !Fddd
1868    symrel(:,:,3) = genmpm(:,:)
1869    symrel(:,:,4) = genpmm(:,:)
1870    symrel(:,:,5) = genmmm(:,:)
1871    symrel(:,:,6) = genppm(:,:)
1872    symrel(:,:,7) = genpmp(:,:)
1873    symrel(:,:,8) = genmpp(:,:)
1874    if (spgorig==1) then
1875      tnons(:,5)=(/0.25d0,0.25d0,0.25d0/)
1876      tnons(:,6)=(/0.25d0,0.25d0,0.25d0/)
1877      tnons(:,7)=(/0.25d0,0.25d0,0.25d0/)
1878      tnons(:,8)=(/0.25d0,0.25d0,0.25d0/)
1879    else
1880      tnons(:,2)=(/0.75d0,0.75d0,0.d0/)
1881      tnons(:,3)=(/0.75d0,0.d0,0.75d0/)
1882      tnons(:,4)=(/0.d0,0.75d0,0.75d0/)
1883 !      JWZ DEBUGGING BEGIN
1884 !      the following lines were present in 5.7 as of Dec 10 2008 but
1885 !      gave wrong results in a spgroup 70 spgorig 2 case (Na2SO4)
1886 !      tnons(:,5)=(/0.25d0,0.25d0,0.d0/) ! original code
1887 !      tnons(:,6)=(/0.25d0,0.d0,0.25d0/) ! original code
1888 !      tnons(:,7)=(/0.d0,0.25d0,0.25d0/) ! original code
1889 !      here are the corrected values of tnons for this case
1890      tnons(:,5)=(/0.0d0,0.0d0,0.0d0/)
1891      tnons(:,6)=(/0.25d0,0.25d0,0.0d0/)
1892      tnons(:,7)=(/0.25d0,0.0d0,0.25d0/)
1893      tnons(:,8)=(/0.0d0,0.25d0,0.25d0/)
1894 !      JWZ DEBUGGING END
1895    end if
1896    if(shubnikov==3)then
1897      if(spgroupma==529)then
1898        symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1
1899      else if(spgroupma==530)then
1900        symafm(3:4)=-1 ; symafm(7:8)=-1
1901      else if(spgroupma==531)then
1902        symafm(5:8)=-1
1903      end if
1904    end if
1905    nogen=0
1906  case (74)                !Imma
1907    tnons(:,2)=(/0.d0,0.5d0,0.d0/)
1908    symrel(:,:,3) = genmpm(:,:)
1909    tnons(:,3)=(/0.d0,0.5d0,0.d0/)
1910    symrel(:,:,4) = genpmm(:,:)
1911  end select
1912 
1913  if (shubnikov==3) then
1914    select case (spgroupma)
1915    case (59,68,80,89,101,113,125,137,146,158,167,174,182,189,197,&
1916 &     205,213,221,226,231,237,243,270,292,296,308,312,324,328,340,344,&
1917 &     358,370,380,384,420,424,444,448,460,464,472,476)
1918      symafm(2)=-1
1919      symafm(4)=-1
1920    case (69,90,102,114,126,147,175,190,198,206,214,244)
1921      symafm(2:3)=-1
1922    case (60,70,81,91,103,115,127,138,148,159,168,176,183,191,199,&
1923 &     207,215,222,227,232,238,245,252,268,269,293,294,309,310,&
1924 &     325,326,341,342,356,357,368,369,381,382,396,397,421,422,436,445,&
1925 &     446,461,462,473,474,484,485,494,495,504,505,524,536,&
1926 &     542,543,551,557,558)
1927      symafm(3:4)=-1
1928    case (251,267,291,295,307,311,323,327,339,343,355,367,379,383,&
1929 &     395,398,419,423,435,443,447,459,463,&
1930 &     471,475,483,486,493,496,503,506,523,535,541,544,550,556,559)
1931      symafm(2:3)=-1
1932 
1933    end select
1934  end if
1935 
1936  if (nogen>1)then
1937    call bldgrp(msym,nogen,nsym,symafm,symrel,tnons)
1938  end if
1939 
1940  call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig)
1941 
1942 !DEBUG
1943 !write(std_out,*)'symsgortho : end of symmetry assignement'
1944 !ENDDEBUG
1945 
1946 end subroutine symsgortho

m_symsg/symsgtetra [ Functions ]

[ Top ] [ m_symsg ] [ Functions ]

NAME

 symsgtetra

FUNCTION

 Yields all the TETRAGONAL symmetry operations starting from the space group symbol.
 according to the International Tables of Crystallography, 1983.

INPUTS

 msym = default number of symmetries
 nsym = the number of symmetry operations
 shubnikov= magnetic type of the space group to be generated
 spgorig = the origin choice (1 or 2) for the axes system
 spgroup = the numeric symbol of the space groups
 spgroupma= number of the magnetic space group

OUTPUT

 symafm(nsym)=(anti)ferromagnetic part of symmetry operations
 symrel(3,3,nsym) = 3D matrix containg symmetry operations
 tnons(3,nsym) = 2D matrix containing translations associated

PARENTS

      gensymspgr

CHILDREN

      bldgrp,spgdata

SOURCE

1978 subroutine symsgtetra(msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma,symafm,symrel,tnons)
1979 
1980 
1981 !This section has been created automatically by the script Abilint (TD).
1982 !Do not modify the following lines by hand.
1983 #undef ABI_FUNC
1984 #define ABI_FUNC 'symsgtetra'
1985 !End of the abilint section
1986 
1987  implicit none
1988 
1989 !Arguments ------------------------------------
1990 !scalars
1991  integer,intent(in) :: msym,nsym,shubnikov,spgaxor,spgorig,spgroup
1992  integer,intent(in) :: spgroupma
1993 !arrays
1994  integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i
1995  real(dp),intent(inout) :: tnons(3,msym) !vz_i
1996 
1997 !Local variables ------------------------------
1998 !integer :: isym
1999 !scalars
2000  integer :: nogen,sporder
2001  character(len=1) :: brvsb
2002  character(len=15) :: intsb,ptintsb,ptschsb,schsb
2003  character(len=35) :: intsbl
2004 !arrays
2005  integer :: gen4m(3,3),genmmm(3,3),genmmp(3,3),genmpm(3,3),genmpp(3,3)
2006  integer :: genpmm(3,3),genpmp(3,3),genppm(3,3)
2007 
2008 ! *************************************************************************
2009 !DEBUG
2010 !write(std_out,*) ' symsgtetra: ',spgroup,shubnikov,spgroupma
2011 !ENDDEBUG
2012  nogen=0
2013 
2014  tnons(:,:)=zero
2015 !The identity operation belongs to all space groups
2016  symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1
2017 
2018 !The next operation belongs to most TETRAGONAL space groups, except 81,82,111,121.
2019  symrel(:,:,2)=0 ; symrel(1,2,2)=-1 ; symrel(2,1,2)=1 ; symrel(3,3,2)=1
2020 
2021 !Predefine some generators
2022  genmpp(:,:)=0 ; genmpp(1,1)=-1 ; genmpp(2,2)= 1 ; genmpp(3,3)= 1
2023  genpmp(:,:)=0 ; genpmp(1,1)= 1 ; genpmp(2,2)=-1 ; genpmp(3,3)= 1
2024  genppm(:,:)=0 ; genppm(1,1)= 1 ; genppm(2,2)= 1 ; genppm(3,3)=-1
2025  genpmm(:,:)=0 ; genpmm(1,1)= 1 ; genpmm(2,2)=-1 ; genpmm(3,3)=-1
2026  genmpm(:,:)=0 ; genmpm(1,1)=-1 ; genmpm(2,2)= 1 ; genmpm(3,3)=-1
2027  genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)= 1
2028  genmmm(:,:)=0 ; genmmm(1,1)=-1 ; genmmm(2,2)=-1 ; genmmm(3,3)=-1
2029  gen4m(:,:)=0  ; gen4m(1,2)=-1 ; gen4m(2,1)=1 ; gen4m(3,3)=-1
2030 
2031 !Default non-magnetic behaviour
2032  symafm(1:nsym)=1
2033 
2034 
2035 !assigns the generators to each space group
2036  select case (spgroup)
2037 !  TETRAGONAL space groups
2038  case (75,79,83,87)        !P4, I4, P4/m, I4/m
2039    nogen=2
2040  case (76,80)                !P41, I41
2041    tnons(:,2)=(/0.d0,0.d0,0.25d0/)
2042    nogen=2
2043  case (77,84)                !P42, P42/m
2044    tnons(:,2)=(/0.d0,0.d0,0.5d0/)
2045    nogen=2
2046  case (78)                !P43
2047    tnons(:,2)=(/0.d0,0.d0,0.75d0/)
2048    nogen=2
2049  case (81,82)                !PB4, IB4
2050    symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1
2051    symrel(:,:,2)=0 ; symrel(1,1,2)=-1; symrel(2,2,2)=-1; symrel(3,3,2)=1
2052    symrel(:,:,3)=0 ; symrel(1,2,3)=1 ; symrel(2,1,3)=-1 ; symrel(3,3,3)=-1
2053    symrel(:,:,4)=0 ; symrel(1,2,4)=-1 ; symrel(2,1,4)=1 ; symrel(3,3,4)=-1
2054    nogen=0
2055    if (shubnikov==3) then
2056      symafm(3:4)=-1
2057    end if
2058  case (85,86,88)                !P4/n
2059    symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 ; symrel(:,:,5) = -1*symrel(:,:,1)
2060    symrel(:,:,2)=0 ; symrel(1,1,2)=-1; symrel(2,2,2)=-1; symrel(3,3,2)=1 ; symrel(:,:,6) = -1*symrel(:,:,2)
2061    symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(3,3,3)=1 ; symrel(:,:,7) = -1*symrel(:,:,3)
2062    symrel(:,:,4)=0 ; symrel(1,2,4)=1 ; symrel(2,1,4)=-1 ; symrel(3,3,4)=1 ; symrel(:,:,8) = -1*symrel(:,:,4)
2063    nogen=0
2064    if(spgorig==1) then
2065      select case (spgroup)
2066      case (85)                !P4/n
2067        tnons(:,3)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
2068        tnons(:,5)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.d0/)
2069      case (86)                !P42/n
2070        tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
2071        tnons(:,5)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.5d0/)
2072      case (88)                 !I41/a
2073        tnons(:,2)=(/0.5d0,0.5d0,0.5d0/)
2074        tnons(:,3)=(/0.0d0,0.5d0,0.25d0/)
2075        tnons(:,4)=(/0.5d0,0.0d0,0.75d0/)
2076        tnons(:,5)=(/0.0d0,0.5d0,0.25d0/)
2077        tnons(:,6)=(/0.5d0,0.0d0,0.75d0/)
2078        tnons(:,7)=(/0.5d0,0.5d0,0.5d0/)
2079      end select
2080    else
2081      select case (spgroup)
2082      case (85)          !P4/n
2083        tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ;  tnons(:,6)=(/0.5d0,0.5d0,0.d0/)
2084        tnons(:,3)=(/0.5d0,0.0d0,0.d0/) ;  tnons(:,7)=(/0.5d0,0.0d0,0.d0/)
2085        tnons(:,4)=(/0.0d0,0.5d0,0.d0/) ;  tnons(:,8)=(/0.0d0,0.5d0,0.d0/)
2086      case (86)          !P42/n
2087        tnons(:,2)=(/0.5d0,0.5d0,0.0d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.0d0/)
2088        tnons(:,3)=(/0.0d0,0.5d0,0.5d0/) ; tnons(:,7)=(/0.0d0,0.5d0,0.5d0/)
2089        tnons(:,4)=(/0.5d0,0.0d0,0.5d0/) ; tnons(:,8)=(/0.5d0,0.0d0,0.5d0/)
2090      case (88)          !I41/a
2091        tnons(:,2)=(/0.5d0,0.0d0,0.5d0/) ;  tnons(:,6)=(/0.5d0,0.0d0,0.5d0/)
2092        tnons(:,3)=(/0.75d0,0.25d0,0.25d0/)
2093        tnons(:,7)=(/0.25d0,0.75d0,0.75d0/)
2094        tnons(:,4)=(/0.75d0,0.75d0,0.75d0/)
2095        tnons(:,8)=(/0.25d0,0.25d0,0.25d0/)
2096      end select
2097    end if
2098    if (shubnikov==3) then
2099      select case (spgroupma)
2100      case(61,69,83)
2101        symafm(3)=-1;symafm(4)=-1;symafm(7)=-1; symafm(8)=-1
2102      case(62,70,84)
2103        symafm(5)=-1;symafm(6)=-1;symafm(7)=-1; symafm(8)=-1
2104      case(63,71,85)
2105        symafm(3)=-1;symafm(5)=-1;symafm(4)=-1; symafm(6)=-1
2106      end select
2107    end if
2108  case (89,97)                !P422, I422
2109    symrel(:,:,3) = genmpm(:,:)
2110    nogen=3
2111  case (90)                !P4212
2112    tnons(:,2)=(/0.5d0,0.5d0,0.0d0/)
2113    symrel(:,:,3) = genmpm(:,:)
2114    tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
2115    nogen=3
2116  case (91)                !P4122
2117    tnons(:,2)=(/0.d0,0.d0,0.25d0/)
2118    symrel(:,:,3) = genmpm(:,:)
2119    nogen=3
2120  case (92)                !P41212
2121    tnons(:,2)=(/0.5d0,0.5d0,0.25d0/)
2122    symrel(:,:,3) = genmpm(:,:)
2123    tnons(:,3)=(/0.5d0,0.5d0,0.25d0/)
2124    nogen=3
2125  case (93)                !P4222
2126    tnons(:,2)=(/0.d0,0.d0,0.5d0/)
2127    symrel(:,:,3) = genmpm(:,:)
2128    nogen=3
2129  case (94)                !P42212
2130    tnons(:,2)=(/0.5d0,0.5d0,0.5d0/)
2131    symrel(:,:,3) = genmpm(:,:)
2132    tnons(:,3)=(/0.5d0,0.5d0,0.5d0/)
2133    nogen=3
2134  case (95)                !P4322
2135    tnons(:,2)=(/0.d0,0.d0,0.75d0/)
2136    symrel(:,:,3) = genmpm(:,:)
2137    nogen=3
2138  case (96)                !P43212
2139    tnons(:,2)=(/0.5d0,0.5d0,0.75d0/)
2140    symrel(:,:,3) = genmpm(:,:)
2141    tnons(:,3)=(/0.5d0,0.5d0,0.75d0/)
2142    nogen=3
2143  case (98)                !I4122
2144    tnons(:,2)=(/0.d0,0.5d0,0.25d0/)
2145    symrel(:,:,3) = genmpm(:,:)
2146    tnons(:,3)=(/0.5d0,0.0d0,0.75d0/)
2147    nogen=3
2148  case (99,107,123,139)        !P4mm, I4mm, P4/mmm, I4/mmm
2149    symrel(:,:,3) = genpmp(:,:)
2150    nogen=3
2151  case (100)                !P4bm
2152    symrel(:,:,3) = genpmp(:,:)
2153    tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
2154    nogen=3
2155  case (101,132)                !P42cm, P42/mcm
2156    tnons(:,2)=(/0.d0,0.d0,0.5d0/)
2157    symrel(:,:,3) = genpmp(:,:)
2158    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
2159    nogen=3
2160  case (102)                !P42nm
2161    tnons(:,2)=(/0.5d0,0.5d0,0.5d0/)
2162    symrel(:,:,3) = genpmp(:,:)
2163    tnons(:,3)=(/0.5d0,0.5d0,0.5d0/)
2164    nogen=3
2165  case (103,124)                !P4cc, P4/mcc
2166    symrel(:,:,3) = genpmp(:,:)
2167    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
2168    nogen=3
2169  case (104)                !P4nc
2170    symrel(:,:,3) = genpmp(:,:)
2171    tnons(:,3)=(/0.5d0,0.5d0,0.5d0/)
2172    nogen=3
2173  case (105,131)                !P42mc, P42/mmc
2174    tnons(:,2)=(/0.d0,0.d0,0.5d0/)
2175    symrel(:,:,3) = genpmp(:,:)
2176    nogen=3
2177  case (106,135)                !P42bc, P42/mbc
2178    tnons(:,2)=(/0.d0,0.d0,0.5d0/)
2179    symrel(:,:,3) = genpmp(:,:)
2180    tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
2181    nogen=3
2182  case (108)                !I4cm
2183    symrel(:,:,3) = genpmp(:,:)
2184    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
2185    nogen=3
2186  case (109)                !I41md
2187    tnons(:,2)=(/0.d0,0.5d0,0.25d0/)
2188    symrel(:,:,3) = genpmp(:,:)
2189    symrel(:,:,4) = genmmp(:,:)
2190    tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
2191    nogen=4
2192  case (110)                !I41cd
2193    tnons(:,2)=(/0.d0,0.5d0,0.25d0/)
2194    symrel(:,:,3) = genpmp(:,:)
2195    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
2196    symrel(:,:,4) = genmmp(:,:)
2197    tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
2198    nogen=4
2199  case (111,121)                !PB42m, IB42m
2200    symrel(:,:,2) = gen4m(:,:)
2201    symrel(:,:,3) = genmpm(:,:)
2202    nogen=3
2203  case (112)                !PB42c
2204    symrel(:,:,2) = gen4m(:,:)
2205    symrel(:,:,3) = genmpm(:,:)
2206    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
2207    nogen=3
2208  case (113)                !PB421m
2209    symrel(:,:,2) = gen4m(:,:)
2210    symrel(:,:,3) = genmpm(:,:)
2211    tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
2212    nogen=3
2213  case (114)                !PB421c
2214    symrel(:,:,2) = gen4m(:,:)
2215    symrel(:,:,3) = genmpm(:,:)
2216    tnons(:,3)=(/0.5d0,0.5d0,0.5d0/)
2217    nogen=3
2218  case (115,119)                !PB4m2,IB4m2
2219    symrel(:,:,2) = gen4m(:,:)
2220    symrel(:,:,3) = genpmp(:,:)
2221    nogen=3
2222  case (116,120)                !PB4c2, IB4c2
2223    symrel(:,:,2) = gen4m(:,:)
2224    symrel(:,:,3) = genpmp(:,:)
2225    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
2226    nogen=3
2227  case (117)                !PB4b2
2228    symrel(:,:,2) = gen4m(:,:)
2229    symrel(:,:,3) = genpmp(:,:)
2230    tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
2231    nogen=3
2232  case (118)                !PB4n2
2233    symrel(:,:,2) = gen4m(:,:)
2234    symrel(:,:,3) = genpmp(:,:)
2235    tnons(:,3)=(/0.5d0,0.5d0,0.5d0/)
2236    nogen=3
2237  case (122)                !IB42d
2238    symrel(:,:,2)=0 ; symrel(1,2,2)=-1 ; symrel(2,1,2)=1 ; symrel(3,3,2)=-1
2239    symrel(:,:,3) = genmpm(:,:)
2240    tnons(:,3)=(/0.5d0,0.d0,0.75d0/)
2241    nogen=3
2242  case (125,126,129,130,133,134,137,138,141,142)
2243    symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1
2244    symrel(:,:,2)=0 ; symrel(1,1,2)=-1; symrel(2,2,2)=-1; symrel(3,3,2)=1
2245    symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(3,3,3)=1
2246    symrel(:,:,4)=0 ; symrel(1,2,4)=1 ; symrel(2,1,4)=-1 ; symrel(3,3,4)=1
2247    symrel(:,:,5) = genmpm(:,:)
2248    symrel(:,:,6) = genmmm(:,:)
2249    if (spgorig==1) then
2250      select case (spgroup)
2251      case (125)                !P4/nbm
2252        tnons(:,6)=(/0.5d0,0.5d0,0.d0/)
2253      case (126)         !P4/nnc
2254        tnons(:,6)=(/0.5d0,0.5d0,0.5d0/)
2255      case (129)                !P4/nmm
2256        tnons(:,3)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
2257        tnons(:,5)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.d0/)
2258      case (130)                !P4/ncc
2259        tnons(:,3)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.d0/)
2260        tnons(:,5)=(/0.5d0,0.5d0,0.5d0/)
2261        tnons(:,6)=(/0.5d0,0.5d0,0.d0/)
2262      case (133)         !P42/nbc
2263        tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
2264        tnons(:,5)=(/0.d0,0.d0,0.5d0/)
2265        tnons(:,6)=(/0.5d0,0.5d0,0.5d0/)
2266      case (134)         !P42/nnm
2267        tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
2268        tnons(:,6)=(/0.5d0,0.5d0,0.5d0/)
2269      case (137)         !P42/nmc
2270        tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
2271        tnons(:,5)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.5d0/)
2272      case (138)         !P42/ncm
2273        tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.5d0/)
2274        tnons(:,5)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.5d0/)
2275      case (141)         !I41/amd
2276        tnons(:,2)=(/0.5d0,0.5d0,0.5d0/)
2277        tnons(:,3)=(/0.d0,0.5d0,0.25d0/)
2278        tnons(:,4)=(/0.5d0,0.d0,0.75d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.75d0/)
2279        tnons(:,6)=(/0.d0,0.5d0,0.25d0/)
2280      case (142)         !I41/acd
2281        tnons(:,2)=(/0.5d0,0.5d0,0.5d0/)
2282        tnons(:,3)=(/0.d0,0.5d0,0.25d0/)
2283        tnons(:,4)=(/0.5d0,0.d0,0.75d0/)
2284        tnons(:,5)=(/0.5d0,0.d0,0.25d0/)
2285        tnons(:,6)=(/0.d0,0.5d0,0.25d0/)
2286      end select
2287    else
2288      select case (spgroup)
2289      case (125)         !P4/nbm
2290        tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.d0/)
2291        tnons(:,4)=(/0.0d0,0.5d0,0.d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.d0/)
2292      case (126)                !P4/nnc
2293        tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.d0/)
2294        tnons(:,4)=(/0.0d0,0.5d0,0.d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.5d0/)
2295      case (129)                !P4/nmm
2296        tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.d0/)
2297        tnons(:,4)=(/0.0d0,0.5d0,0.d0/) ;  tnons(:,5)=(/0.0d0,0.5d0,0.d0/)
2298      case (130)         !P4/ncc
2299        tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ;  tnons(:,3)=(/0.5d0,0.d0,0.d0/)
2300        tnons(:,4)=(/0.0d0,0.5d0,0.d0/) ; tnons(:,5)=(/0.0d0,0.5d0,0.5d0/)
2301      case (133)        !P42/nbc
2302        tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
2303        tnons(:,4)=(/0.0d0,0.5d0,0.5d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.d0/)
2304      case (134)                !P42/nnm
2305        tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
2306        tnons(:,4)=(/0.0d0,0.5d0,0.5d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.5d0/)
2307      case (137)         !P42/nmc
2308        tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
2309        tnons(:,4)=(/0.0d0,0.5d0,0.5d0/) ; tnons(:,5)=(/0.d0,0.5d0,0.d0/)
2310      case (138)         !P42/ncm
2311        tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.5d0/)
2312        tnons(:,4)=(/0.0d0,0.5d0,0.5d0/) ; tnons(:,5)=(/0.d0,0.5d0,0.5d0/)
2313      case (141)         !I41/amd
2314        tnons(:,2)=(/0.5d0,0.d0,0.5d0/) ; tnons(:,3)=(/0.25d0,0.75d0,0.25d0/)
2315        tnons(:,4)=(/0.25d0,0.25d0,0.75d0/) ;  tnons(:,5)=(/0.5d0,0.d0,0.5d0/)
2316      case (142)         !I41/acd
2317        tnons(:,2)=(/0.5d0,0.d0,0.5d0/) ; tnons(:,3)=(/0.25d0,0.75d0,0.25d0/)
2318        tnons(:,4)=(/0.25d0,0.25d0,0.75d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.d0/)
2319      end select
2320    end if
2321    nogen=6
2322  case (127)                !P4/mbm
2323    symrel(:,:,3) = genmpm(:,:)
2324    tnons(:,3)=(/0.5d0,0.5d0,0.d0/)
2325    nogen=3
2326  case (128)                !P4/mnc
2327    symrel(:,:,3) = genmpm(:,:)
2328    tnons(:,3)=(/0.5d0,0.5d0,0.5d0/)
2329    nogen=3
2330  case (136)                !P42/mnm
2331    tnons(:,2)=(/0.5d0,0.5d0,0.5d0/)
2332    symrel(:,:,3) = genmpm(:,:)
2333    tnons(:,3)=(/0.5d0,0.5d0,0.5d0/)
2334    nogen=3
2335  case (140)                !I4/mcm
2336    symrel(:,:,3) = genmpm(:,:)
2337    tnons(:,3)=(/0.d0,0.d0,0.5d0/)
2338    nogen=3
2339  end select
2340 
2341  if(shubnikov==3)then
2342    select case(spgroupma)
2343    case(3,9,15,21,27,31,45,47,53,55,77,79,&
2344 &     89,97,105,113,121,129,137,145,153,159,166,174,182,190,198,206,&
2345 &     214,222,230,236,242,248,254,262,270,278,285,293,301,309,317,&
2346 &     323,330,336,343,346,355,358,391,392,403,404,&
2347 &     439,442,451,454,487,490,499,500,535,&
2348 &     538,545,546)
2349      symafm(2)=-1
2350    case(90,98,106,114,122,130,138,146,154,160,167,175,183,191,199,&
2351 &     207,215,223,231,237,243,249,255,263,271,279,287,295,303,311,319,&
2352 &     325,331,337,345,347,357,359,389,393,401,405,441,443,453,455,&
2353 &     489,491,497,501,537,539,543,547)
2354      symafm(3)=-1
2355    case(91,99,107,115,123,131,139,147,155,161,165,173,181,189,&
2356 &     197,205,213,221,229,235,241,247,253,261,269,277,286,294,302,310,&
2357 &     318,324,329,335,342,344,354,356,390,394,402,406,&
2358 &     438,440,450,452,486,488,498,502,534,536,544,548)
2359      symafm(2:3)=-1
2360    case(365,377,413,425,461,473,509,521)
2361      symafm(2)=-1
2362      symafm(5:6)=-1
2363    case(366,378,414,426,462,474,510,522,554,564)
2364      symafm(3:5)=-1
2365    case(367,379,415,427,463,475,511,523,555,565)
2366      symafm(3:4)=-1
2367    case(368,380,416,428,464,476,512,524,556,566)
2368      symafm(3:4)=-1
2369      symafm(6)=-1
2370    case(369,381,417,429,465,477,513,525)
2371      symafm(2)=-1
2372      symafm(5)=-1
2373    case(370,382,418,430,466,478,514,526,558,568)
2374      symafm(3:6)=-1
2375    case(371,383,419,431,467,479,515,527,559,569)
2376      symafm(6)=-1
2377    case(553,563)
2378      symafm(5:6)=-1
2379    case(557,567)
2380      symafm(5)=-1
2381    end select
2382  end if
2383 
2384 !DEBUG
2385 !write(std_out,*)' symsgtetra : symrel(:,:,6)=',symrel(:,:,6)
2386 !ENDDEBUG
2387 
2388  if (nogen>1) then
2389    call bldgrp(msym,nogen,nsym,symafm,symrel,tnons)
2390  end if
2391 
2392  call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig)
2393 
2394 !DEBUG
2395 !write(std_out,*)' symsgtetra : exit'
2396 !do isym=1,nsym
2397 !write(std_out,'(i3,2x,9i3,3es13.3,i3)') isym,symrel(:,:,isym),tnons(:,isym),symafm(isym)
2398 !end do
2399 !ENDDEBUG
2400 
2401 end subroutine symsgtetra