TABLE OF CONTENTS
- ABINIT/m_symsg
- m_symsg/bldgrp
- m_symsg/symsgcube
- m_symsg/symsghexa
- m_symsg/symsgmono
- m_symsg/symsgortho
- m_symsg/symsgtetra
ABINIT/m_symsg [ Modules ]
NAME
m_symsg
FUNCTION
COPYRIGHT
Copyright (C) 1999-2024 ABINIT group (RC,XG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_symsg 23 24 use defs_basis 25 use m_errors 26 use m_abicore 27 28 use m_time, only : timab 29 use m_spgdata, only : spgdata 30 31 implicit none 32 33 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
SOURCE
2362 subroutine bldgrp(msym,nogen,nsym,symafm,symrel,tnons) 2363 2364 implicit none 2365 2366 !Arguments ------------------------------------ 2367 !scalars 2368 integer,intent(in) :: msym,nsym 2369 integer,intent(inout) :: nogen 2370 !arrays 2371 integer,intent(inout) :: symafm(msym),symrel(3,3,msym) 2372 real(dp),intent(inout) :: tnons(3,msym) 2373 2374 !Local variables ------------------------------ 2375 !matrintoper(3,3) & matrinttransl(3) are intermediate arrays of the new 2376 ! symmetry operations obtained, in order to check their uniqueness. 2377 !flagop,flagtr = flags used during the checking of the similarity between 2378 ! the obtained operation and the already existent ones 2379 !ii,ijk,ijkl,jjj,kk = counters in the cycles 2380 !scalars 2381 integer :: flagop,flagtr,ii,ijk,ijkl,jjj,kk,matrintsymafm,nogen_new 2382 real(dp) :: nastyzero 2383 character(len=500) :: message 2384 !arrays 2385 integer :: bcksymafm(2*msym),bcksymrel(3,3,2*msym),matrintoper(3,3) 2386 real(dp) :: bcktnons(3,2*msym),matrinttransl(3) 2387 2388 ! ************************************************************************* 2389 2390 nastyzero=0.1 2391 2392 !DEBUG 2393 write(std_out,*)' bldgrp : enter, builds the space group symmetry ' 2394 write(std_out,*)' bldgrp : number of generators : ',nogen 2395 write(std_out,*)' bldgrp : nsym,msym=',nsym,msym 2396 !ENDDEBUG 2397 2398 if (nogen<1) then 2399 write(message, '(a,i4,a,a,a,a,a)' )& 2400 & 'The number of generators nogen is ',nogen,& 2401 & 'and it should be greater than one',ch10,& 2402 & 'This is not allowed. ',ch10,& 2403 & 'Action: Contact ABINIT group ' 2404 ABI_ERROR(message) 2405 end if 2406 2407 !Transfer the generators to bcksymrel 2408 do ii=1,nogen 2409 bcksymrel(:,:,ii)=symrel(:,:,ii) 2410 bcktnons(:,ii)=tnons(:,ii) 2411 bcksymafm(ii)=symafm(ii) 2412 end do 2413 2414 !DEBUG 2415 write(std_out,*)' Describe the different generators (index,symrel,tnons,symafm)' 2416 do ii=1,nogen 2417 write(std_out,'(i3,2x,9i3,3es12.2,i3)')ii,symrel(:,:,ii),tnons(:,ii),symafm(ii) 2418 end do 2419 !ENDDEBUG 2420 2421 !Simply iterate until the group is complete 2422 do ijkl=1,nsym 2423 2424 ! DEBUG 2425 write(std_out,*)' bldgrp : in loop, ijkl,nogen=',ijkl,nogen 2426 ! ENDDEBUG 2427 2428 nogen_new=nogen 2429 2430 do jjj=2,nogen 2431 do kk=2,nogen 2432 2433 ! Computing block of the new symmetry operation according to: 2434 ! ! $ { R1 | v1 }{ R2 | v2 } = { R1.R2 | v1+R1.v2 } $ 2435 matrintoper(:,:) = matmul(bcksymrel(:,:,jjj),bcksymrel(:,:,kk)) 2436 matrinttransl(:) = bcktnons(:,jjj)+matmul(bcksymrel(:,:,jjj),bcktnons(:,kk)) 2437 matrintsymafm = bcksymafm(jjj)*bcksymafm(kk) 2438 2439 ! Rescaling translation between 0 and 1 2440 do ii=1,3 2441 if (matrinttransl(ii)>=0.9) then 2442 do while (matrinttransl(ii)>=0.9) 2443 matrinttransl(ii)=matrinttransl(ii)-1.0 2444 end do 2445 end if 2446 if (matrinttransl(ii)<0.0) then 2447 do while (matrinttransl(ii)<0.0) 2448 matrinttransl(ii)=matrinttransl(ii)+1.0 2449 end do 2450 end if 2451 if ( abs(matrinttransl(ii))<nastyzero) matrinttransl(ii)=0.0 2452 if ( abs(matrinttransl(ii)-1.0)<nastyzero) matrinttransl(ii)=0.0 2453 end do 2454 2455 ! Cheking block to validate the new symmetry operation 2456 do ijk=1,nogen_new 2457 2458 flagop=0 ; flagtr=0 2459 2460 ! Check for rotation similarity 2461 if(sum((matrintoper-bcksymrel(:,:,ijk))**2)==0)flagop=1 2462 2463 ! Check for translation similarity 2464 if(maxval((matrinttransl-bcktnons(:,ijk))**2)<nastyzero**2)flagtr=1 2465 2466 if(flagop+flagtr==2)exit 2467 2468 end do 2469 2470 ! Add the new determined symmetry if it is unique 2471 if (flagtr+flagop<2) then 2472 nogen_new=nogen_new+1 2473 bcksymrel(:,:,nogen_new)=matrintoper(:,:) 2474 bcktnons(:,nogen_new)=matrinttransl(:) 2475 bcksymafm(nogen_new)=matrintsymafm 2476 ! DEBUG 2477 write(std_out,*)' added one more symmetry : nogen_new=',nogen_new 2478 write(std_out,'(i3,2x,9i3,3es12.2,i3)')& 2479 & nogen_new,bcksymrel(:,:,nogen_new),bcktnons(:,nogen_new),bcksymafm(nogen_new) 2480 ! ENDDEBUG 2481 end if 2482 2483 end do 2484 end do 2485 2486 nogen=nogen_new 2487 2488 if(nogen==nsym)exit 2489 2490 end do 2491 2492 !Transfer of the calculated symmetry to the routine output 2493 if (nogen==nsym) then 2494 symrel(:,:,1:nsym)=bcksymrel(:,:,1:nsym) 2495 tnons(:,1:nsym)=bcktnons(:,1:nsym) 2496 symafm(1:nsym)=bcksymafm(1:nsym) 2497 else 2498 ! Problem with the generation of the symmetry operations 2499 write(message, '(a,i7,a,a,i7)' )& 2500 & 'The symmetries obtained are ',nogen,ch10,& 2501 & 'and they should be ',nsym 2502 ABI_BUG(message) 2503 end if 2504 2505 !DEBUG 2506 write(std_out,*)' bldgrp : exit with ',nogen,' operation symmetries' 2507 !ENDDEBUG 2508 2509 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
SOURCE
70 subroutine symsgcube(msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma,symafm,symrel,tnons) 71 72 implicit none 73 74 !Arguments ------------------------------------ 75 !scalars 76 integer,intent(in) :: msym,shubnikov,spgaxor,spgorig,spgroup,spgroupma 77 integer,intent(inout) :: nsym !vz_i 78 !arrays 79 integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i 80 real(dp),intent(out) :: tnons(3,msym) 81 82 !Local variables ----------------------------- 83 !scalars 84 integer :: ii,nogen,sporder 85 character(len=1) :: brvsb 86 character(len=15) :: intsb,ptintsb,ptschsb,schsb 87 character(len=35) :: intsbl 88 !arrays 89 integer :: gen1(3,3),gen2(3,3),gen3(3,3),gen4(3,3),gen5(3,3),gen6(3,3) 90 integer :: gen7(3,3),gen8(3,3),gen9(3,3),genmmm(3,3),genmmp(3,3),genmpm(3,3) 91 integer :: genmpp(3,3),genpmm(3,3),genpmp(3,3),genppm(3,3),genrot(3,3) 92 integer :: genswm(3,3),genswp(3,3) 93 real(dp) :: tsec(2) 94 95 !************************************************************************* 96 97 !DEBUG 98 !write(std_out,*) ' symsgcube : enter with spgroup ',spgroup,' and ',' origin choice ',spgorig 99 !write(std_out,*) ' msym,nsym=',msym,nsym 100 !ENDDEBUG 101 102 103 !The identity operation belongs to all space groups 104 symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 105 106 !Initialize the associated translations matrix to 0 107 do ii=1,msym 108 tnons(:,ii)= 0.0d0 109 end do 110 nogen=0 111 112 !Predefine some generators 113 genmpp(:,:)=0 ; genmpp(1,1)=-1 ; genmpp(2,2)= 1 ; genmpp(3,3)= 1 114 genpmp(:,:)=0 ; genpmp(1,1)= 1 ; genpmp(2,2)=-1 ; genpmp(3,3)= 1 115 genppm(:,:)=0 ; genppm(1,1)= 1 ; genppm(2,2)= 1 ; genppm(3,3)=-1 116 genpmm(:,:)=0 ; genpmm(1,1)= 1 ; genpmm(2,2)=-1 ; genpmm(3,3)=-1 117 genmpm(:,:)=0 ; genmpm(1,1)=-1 ; genmpm(2,2)= 1 ; genmpm(3,3)=-1 118 genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)= 1 119 genmmm(:,:)=0 ; genmmm(1,1)=-1 ; genmmm(2,2)=-1 ; genmmm(3,3)=-1 120 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/) ) 121 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/) ) 122 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/) ) 123 124 !Because of problems with the IBM compiler, that does not like reshape 125 !operations, define 9 basic matrices 126 gen1(:,:)=0 ; gen1(1,1)=1 127 gen2(:,:)=0 ; gen2(1,2)=1 128 gen3(:,:)=0 ; gen3(1,3)=1 129 gen4(:,:)=0 ; gen4(2,1)=1 130 gen5(:,:)=0 ; gen5(2,2)=1 131 gen6(:,:)=0 ; gen6(2,3)=1 132 gen7(:,:)=0 ; gen7(3,1)=1 133 gen8(:,:)=0 ; gen8(3,2)=1 134 gen9(:,:)=0 ; gen9(3,3)=1 135 136 !Default non-magnetic behaviour 137 symafm(1:msym)=1 138 139 !************************************************************************* 140 141 !Treat CUBIC groups 142 143 if(195<=spgroup .and. spgroup<=230)then 144 145 select case(spgroup) 146 case (195,196,197) !P23, F23, I23 147 symrel(:,:,2) = genrot(:,:) 148 symrel(:,:,3) = genmmp(:,:) 149 symrel(:,:,4) = genmpm(:,:) 150 symrel(:,:,5) = genpmm(:,:) 151 nogen=5 152 case (200,202,204) !PmB3, FmB3, ImB3 153 symrel(:,:,2) = genrot(:,:) 154 symrel(:,:,3) = genmmp(:,:) 155 symrel(:,:,4) = genmpm(:,:) 156 nogen=4 157 case (198,199) !P213, I213 158 symrel(:,:,2) = genrot(:,:) 159 symrel(:,:,3) = genmmp(:,:) 160 tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 161 symrel(:,:,4) = genmpm(:,:) 162 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 163 symrel(:,:,5) = genpmm(:,:) 164 tnons(:,5)=(/0.5d0,0.5d0,0.d0/) 165 symrel(:,:,6) = gen3(:,:)-gen4(:,:)-gen8(:,:) 166 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 167 symrel(:,:,7) = -gen3(:,:)-gen4(:,:)+gen8(:,:) 168 tnons(:,7)=(/0.5d0,0.d0,0.5d0/) 169 symrel(:,:,8) = -gen3(:,:)+gen4(:,:)-gen8(:,:) 170 tnons(:,8)=(/0.d0,0.5d0,0.5d0/) 171 symrel(:,:,10) = -gen2(:,:)+gen6(:,:)-gen7(:,:) 172 tnons(:,10)=(/0.d0,0.5d0,0.5d0/) 173 symrel(:,:,11) = gen2(:,:)-gen6(:,:)-gen7(:,:) 174 tnons(:,11)=(/0.5d0,0.5d0,0.0d0/) 175 symrel(:,:,12) = -gen2(:,:)-gen6(:,:)+gen7(:,:) 176 tnons(:,12)=(/0.5d0,0.d0,0.5d0/) 177 symrel(:,:,9) = gen2(:,:)+gen6(:,:)+gen7(:,:) 178 case (201) !Pn-3 179 if (spgorig==1) then 180 symrel(:,:,2) = genrot(:,:) 181 symrel(:,:,3) = genmmp(:,:) 182 symrel(:,:,4) = genmpm(:,:) 183 symrel(:,:,5) = genmmm(:,:) 184 tnons(:,5)=(/0.5d0,0.5d0,0.5d0/) 185 nogen=5 186 if(shubnikov==3)symafm(5)=-1 187 call bldgrp(msym,nogen,nsym,symafm,symrel,tnons) 188 else 189 symrel(:,:,2) = genrot(:,:) 190 symrel(:,:,3) = genmmp(:,:) 191 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 192 symrel(:,:,4) = genmpm(:,:) 193 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 194 symrel(:,:,5) = genpmm(:,:) 195 tnons(:,5)=(/0.d0,0.5d0,0.5d0/) 196 symrel(:,:,6) = gen3(:,:)-gen4(:,:)-gen8(:,:) 197 tnons(:,6)=(/0.d0,0.5d0,0.5d0/) 198 symrel(:,:,7) = -gen3(:,:)-gen4(:,:)+gen8(:,:) 199 tnons(:,7)=(/0.5d0,0.5d0,0.d0/) 200 symrel(:,:,8) = -gen3(:,:)+gen4(:,:)-gen8(:,:) 201 tnons(:,8)=(/0.5d0,0.d0,0.5d0/) 202 symrel(:,:,10) = -gen2(:,:)+gen6(:,:)-gen7(:,:) 203 tnons(:,10)=(/0.5d0,0.d0,0.5d0/) 204 symrel(:,:,11) = gen2(:,:)-gen6(:,:)-gen7(:,:) 205 tnons(:,11)=(/0.d0,0.5d0,0.5d0/) 206 symrel(:,:,12) = -gen2(:,:)-gen6(:,:)+gen7(:,:) 207 tnons(:,12)=(/0.5d0,0.5d0,0.d0/) 208 symrel(:,:,9) = gen2(:,:)+gen6(:,:)+gen7(:,:) 209 do ii=1,12 210 symrel(:,:,ii+12)= - symrel(:,:,ii) 211 tnons(:,ii+12)=tnons(:,ii) 212 if(shubnikov==3)symafm(ii+12)=-1 213 end do 214 end if 215 nogen=0 216 case (205,206) !Pa-3,Ia-3 217 symrel(:,:,2) = genrot(:,:) 218 symrel(:,:,3) = genmmp(:,:) 219 tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 220 symrel(:,:,4) = genmpm(:,:) 221 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 222 symrel(:,:,5) = genpmm(:,:) 223 tnons(:,5)=(/0.5d0,0.5d0,0.d0/) 224 symrel(:,:,6) = gen3(:,:)-gen4(:,:)-gen8(:,:) 225 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 226 symrel(:,:,7) = -gen3(:,:)-gen4(:,:)+gen8(:,:) 227 tnons(:,7)=(/0.5d0,0.d0,0.5d0/) 228 symrel(:,:,8) = -gen3(:,:)+gen4(:,:)-gen8(:,:) 229 tnons(:,8)=(/0.0d0,0.5d0,0.5d0/) 230 symrel(:,:,9) = gen2(:,:)+gen6(:,:)+gen7(:,:) 231 symrel(:,:,10) = -gen2(:,:)+gen6(:,:)-gen7(:,:) 232 tnons(:,10)=(/0.d0,0.5d0,0.5d0/) 233 symrel(:,:,11) = gen2(:,:)-gen6(:,:)-gen7(:,:) 234 tnons(:,11)=(/0.5d0,0.5d0,0.d0/) 235 symrel(:,:,12) = -gen2(:,:)-gen6(:,:)+gen7(:,:) 236 tnons(:,12)=(/0.5d0,0.d0,0.5d0/) 237 nogen=0 238 case (203) !FdB3 239 if (spgorig==1) then 240 symrel(:,:,2) = genrot(:,:) 241 symrel(:,:,3) = genmmp(:,:) 242 nogen=3 243 nsym=12 244 call bldgrp(msym,nogen,nsym,symafm,symrel,tnons) 245 do ii=1,12 246 symrel(:,:,ii+12) = - symrel(:,:,ii) 247 tnons(:,ii+12)=(/0.25d0,0.25d0,0.25d0/) 248 if(shubnikov==3)symafm(ii+12)=-1 249 end do 250 nsym=24 251 nogen=0 252 else 253 symrel(:,:,2) = genrot(:,:) 254 symrel(:,:,3) = genmmp(:,:) 255 tnons(:,3)=(/0.25d0,0.25d0,0.d0/) 256 symrel(:,:,4) = genmpm(:,:) 257 tnons(:,4)=(/0.25d0,0.d0,0.25d0/) 258 symrel(:,:,5) = genpmm(:,:) 259 tnons(:,5)=(/0.d0,0.25d0,0.25d0/) 260 symrel(:,:,6) = gen3(:,:)-gen4(:,:)-gen8(:,:) 261 tnons(:,6)=(/0.d0,0.25d0,0.25d0/) 262 symrel(:,:,7) = -gen3(:,:)-gen4(:,:)+gen8(:,:) 263 tnons(:,7)=(/0.25d0,0.25d0,0.d0/) 264 symrel(:,:,8) = -gen3(:,:)+gen4(:,:)-gen8(:,:) 265 tnons(:,8)=(/0.25d0,0.d0,0.25d0/) 266 symrel(:,:,10) = -gen2(:,:)+gen6(:,:)-gen7(:,:) 267 tnons(:,10)=(/0.25d0,0.d0,0.25d0/) 268 symrel(:,:,11) = gen2(:,:)-gen6(:,:)-gen7(:,:) 269 tnons(:,11)=(/0.d0,0.25d0,0.25d0/) 270 symrel(:,:,12) = -gen2(:,:)-gen6(:,:)+gen7(:,:) 271 tnons(:,12)=(/0.25d0,0.25d0,0.d0/) 272 symrel(:,:,9) = gen2(:,:)+gen6(:,:)+gen7(:,:) 273 do ii=1,12 274 symrel(:,:,ii+12) = - symrel(:,:,ii) 275 tnons(:,ii+12)=tnons(:,ii) 276 if(shubnikov==3)symafm(ii+12)=-1 277 end do 278 end if 279 nsym=24 280 nogen=0 281 case (207,209,211) !P432, F432, I432, PmB3m, FmB3m, ImB3m 282 symrel(:,:,2) = genrot(:,:) 283 symrel(:,:,3) = genmmp(:,:) 284 symrel(:,:,4) = genmpm(:,:) 285 symrel(:,:,5) = genswm(:,:) 286 if(shubnikov==3)symafm(5)=-1 287 nogen=5 288 case (208) !P4232 289 symrel(:,:,2) = genrot(:,:) 290 symrel(:,:,3) = genmmp(:,:) 291 symrel(:,:,4) = genmpm(:,:) 292 symrel(:,:,5) = genswm(:,:) 293 tnons(:,5)=(/0.5d0,0.5d0,0.5d0/) 294 if(shubnikov==3)symafm(5)=-1 295 nogen=5 296 case (210) !F4132 297 symrel(:,:,2) = genrot(:,:) 298 symrel(:,:,3) = genmmp(:,:) 299 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 300 symrel(:,:,4) = genmpm(:,:) 301 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 302 symrel(:,:,5) = genswm(:,:) 303 tnons(:,5)=(/0.75d0,0.25d0,0.75d0/) 304 if(shubnikov==3)symafm(5)=-1 305 nogen=5 306 case (212) !P4332 307 symrel(:,:,2) = genrot(:,:) 308 symrel(:,:,3) = genmmp(:,:) 309 tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 310 symrel(:,:,4) = genmpm(:,:) 311 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 312 symrel(:,:,5) = genswm(:,:) 313 tnons(:,5)=(/0.25d0,0.75d0,0.75d0/) 314 if(shubnikov==3)symafm(5)=-1 315 nogen=5 316 case (213,214) !P4132, I4132 317 symrel(:,:,2) = genrot(:,:) 318 symrel(:,:,3) = genmmp(:,:) 319 tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 320 symrel(:,:,4) = genmpm(:,:) 321 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 322 symrel(:,:,5) = genswm(:,:) 323 tnons(:,5)=(/0.75d0,0.25d0,0.25d0/) 324 if(shubnikov==3)symafm(5)=-1 325 nogen=5 326 case (215,216,217) !PB43m, FB43m, IB43m 327 symrel(:,:,2) = genrot(:,:) 328 symrel(:,:,3) = genmmp(:,:) 329 symrel(:,:,4) = genmpm(:,:) 330 symrel(:,:,5) = genswp(:,:) 331 if(shubnikov==3)symafm(5)=-1 332 nogen=5 333 case (218,219) !PB43n, FB343c 334 symrel(:,:,2) = genrot(:,:) 335 symrel(:,:,3) = genmmp(:,:) 336 symrel(:,:,4) = genmpm(:,:) 337 symrel(:,:,5) = genswp(:,:) 338 tnons(:,5)=(/0.5d0,0.5d0,0.5d0/) 339 if(shubnikov==3)symafm(5)=-1 340 nogen=5 341 case (220) !IB43d 342 symrel(:,:,2) = genrot(:,:) 343 symrel(:,:,3) = genmmp(:,:) 344 tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 345 symrel(:,:,4) = genmpm(:,:) 346 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 347 symrel(:,:,5) = genswp(:,:) 348 tnons(:,5)=(/0.25d0,0.25d0,0.25d0/) 349 if(shubnikov==3)symafm(5)=-1 350 nogen=5 351 case (221,225,229) !Pm3m,Fm3m,Im3m 352 symrel(:,:,2) = genrot(:,:) 353 symrel(:,:,3) = genmmp(:,:) 354 symrel(:,:,4) = genmpm(:,:) 355 symrel(:,:,5) = genswm(:,:) 356 if(shubnikov==3)then 357 if(spgroupma==94 .or. spgroupma==95 .or. spgroupma==118 .or. & 358 & spgroupma==119 .or. spgroupma==142 .or. spgroupma==143 )symafm(5)=-1 359 end if 360 nogen=5 361 case (222) !PnB3n 362 if (spgorig==1) then 363 symrel(:,:,2) = genrot(:,:) 364 symrel(:,:,3) = genmmp(:,:) 365 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 366 symrel(:,:,4) = genmpm(:,:) 367 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 368 symrel(:,:,5) = genswm(:,:) 369 tnons(:,5)=(/0.d0,0.d0,0.5d0/) 370 if(shubnikov==3)then 371 if(spgroupma==100 .or. spgroupma==101)symafm(5)=-1 372 end if 373 nogen=5 374 nsym=24 375 call bldgrp(msym,nogen,nsym,symafm,symrel,tnons) 376 do ii=1,24 377 symrel(:,:,ii+24)=-symrel(:,:,ii) 378 tnons(:,ii+24)=tnons(:,ii) 379 if(shubnikov==3)then 380 if(spgroupma==100 .or. spgroupma==102)symafm(ii+24)=-symafm(ii) 381 if(spgroupma==101)symafm(ii+24)=symafm(ii) 382 end if 383 end do 384 nogen=0 ; nsym=48 385 else 386 symrel(:,:,2) = genrot(:,:) 387 symrel(:,:,3) = genmmp(:,:) 388 symrel(:,:,4) = genmpm(:,:) 389 symrel(:,:,5) = genswm(:,:) 390 if(shubnikov==3)then 391 if(spgroupma==100 .or. spgroupma==101)symafm(5)=-1 392 end if 393 nogen=5 394 nsym=24 395 call bldgrp(msym,nogen,nsym,symafm,symrel,tnons) 396 do ii=1,24 397 symrel(:,:,ii+24)=-symrel(:,:,ii) 398 tnons(:,ii+24)=(/0.5d0,0.5d0,0.5d0/) 399 if(shubnikov==3)then 400 if(spgroupma==100 .or. spgroupma==102)symafm(ii+24)=-symafm(ii) 401 if(spgroupma==101)symafm(ii+24)=symafm(ii) 402 end if 403 end do 404 nogen=0 ; nsym=48 405 end if 406 case (223,226) ! PmB3n, FmB3c 407 symrel(:,:,2) = genrot(:,:) 408 symrel(:,:,3) = genmmp(:,:) 409 symrel(:,:,4) = genmpm(:,:) 410 symrel(:,:,5) = genswm(:,:) 411 tnons(:,5)=(/0.5d0,0.5d0,0.5d0/) 412 if(shubnikov==3)then 413 if(spgroupma==106 .or. spgroupma==124 .or. & 414 & spgroupma==107 .or. spgroupma==125 )symafm(5)=-1 415 end if 416 nogen=5 417 case (224) !PnB3m 418 if (spgorig==1) then 419 symrel(:,:,2) = genrot(:,:) 420 symrel(:,:,3) = genmmp(:,:) 421 symrel(:,:,4) = genmpm(:,:) 422 symrel(:,:,5) = genswm(:,:) 423 tnons(:,5)=(/0.5d0,0.5d0,0.5d0/) 424 symrel(:,:,6) = genmmm(:,:) 425 tnons(:,6)=(/0.5d0,0.5d0,0.5d0/) 426 else 427 symrel(:,:,2) = genrot(:,:) 428 symrel(:,:,3) = genmmp(:,:) 429 tnons(:,3)=(/0.5d0,0.5d0,0.0d0/) 430 symrel(:,:,4) = genmpm(:,:) 431 tnons(:,4)=(/0.5d0,0.0d0,0.5d0/) 432 symrel(:,:,5) = genswm(:,:) 433 tnons(:,5)=(/0.5d0,0.5d0,0.d0/) 434 symrel(:,:,6) = genmmm(:,:) 435 end if 436 if(shubnikov==3)then 437 if(spgroupma==112 .or. spgroupma==113)symafm(5)=-1 438 if(spgroupma==112 .or. spgroupma==114)symafm(6)=-1 439 end if 440 nogen=6 441 case (227) !FdB3m 442 if (spgorig==1) then 443 symrel(:,:,2) = genrot(:,:) 444 symrel(:,:,3) = genmmp(:,:) 445 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 446 symrel(:,:,4) = genmpm(:,:) 447 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 448 symrel(:,:,5) = genswm(:,:) 449 tnons(:,5)=(/0.75d0,0.25d0,0.75d0/) 450 symrel(:,:,6) = genmmm(:,:) 451 tnons(:,6)=(/0.25d0,0.25d0,0.25d0/) 452 else 453 symrel(:,:,2) = genrot(:,:) 454 symrel(:,:,3) = genmmp(:,:) 455 tnons(:,3)=(/0.75d0,0.25d0,0.5d0/) 456 symrel(:,:,4) = genmpm(:,:) 457 tnons(:,4)=(/0.25d0,0.5d0,0.75d0/) 458 symrel(:,:,5) = genswm(:,:) 459 tnons(:,5)=(/0.75d0,0.25d0,0.5d0/) 460 symrel(:,:,6) = genmmm(:,:) 461 end if 462 if(shubnikov==3)then 463 if(spgroupma==130 .or. spgroupma==131)symafm(5)=-1 464 if(spgroupma==130 .or. spgroupma==132)symafm(6)=-1 465 end if 466 nogen=6 467 case (228) !FdB3c 468 if (spgorig==1) then 469 symrel(:,:,2) = genrot(:,:) 470 symrel(:,:,3) = genmmp(:,:) 471 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 472 symrel(:,:,4) = genmpm(:,:) 473 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 474 symrel(:,:,5) = genswm(:,:) 475 tnons(:,5)=(/0.75d0,0.25d0,0.75d0/) 476 symrel(:,:,6) = genmmm(:,:) 477 tnons(:,6)=(/0.75d0,0.75d0,0.75d0/) 478 else 479 symrel(:,:,2) = genrot(:,:) 480 symrel(:,:,3) = genmmp(:,:) 481 tnons(:,3)=(/0.25d0,0.75d0,0.5d0/) 482 symrel(:,:,4) = genmpm(:,:) 483 tnons(:,4)=(/0.75d0,0.5d0,0.25d0/) 484 symrel(:,:,5) = genswm(:,:) 485 tnons(:,5)=(/0.75d0,0.25d0,0.d0/) 486 symrel(:,:,6) = genmmm(:,:) 487 end if 488 if(shubnikov==3)then 489 if(spgroupma==136 .or. spgroupma==137)symafm(5)=-1 490 if(spgroupma==136 .or. spgroupma==138)symafm(6)=-1 491 end if 492 nogen=6 493 case (230) !IaB3d 494 symrel(:,:,2) = genrot(:,:) 495 symrel(:,:,3) = genmmp(:,:) 496 tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 497 symrel(:,:,4) = genmpm(:,:) 498 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 499 symrel(:,:,5) = genswm(:,:) 500 tnons(:,5)=(/0.75d0,0.25d0,0.25d0/) 501 if(shubnikov==3)then 502 if(spgroupma==147 .or. spgroupma==148)symafm(5)=-1 503 end if 504 nogen=5 505 end select 506 507 ! End CUBIC space groups 508 end if 509 510 !*************************************************************************** 511 512 call timab(47,1,tsec) 513 514 if (nogen>0) then 515 call bldgrp(msym,nogen,nsym,symafm,symrel,tnons) 516 end if 517 518 call timab(47,2,tsec) 519 520 call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig) 521 522 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
SOURCE
551 subroutine symsghexa(brvltt,msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma,symafm,symrel,tnons) 552 553 implicit none 554 555 !Arguments ------------------------------------ 556 !scalars 557 integer,intent(in) :: msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma 558 integer,intent(inout) :: brvltt !vz_i 559 !arrays 560 integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i 561 real(dp),intent(inout) :: tnons(3,msym) !vz_i 562 563 !Local variables ----------------------------- 564 !scalars 565 integer :: ii,nogen,sporder 566 real(dp),parameter :: fivesixth=5.0d0/6.0d0,twothird=2.0d0/3.0d0 567 character(len=1) :: brvsb 568 character(len=15) :: intsb,ptintsb,ptschsb,schsb 569 character(len=35) :: intsbl 570 !arrays 571 integer :: genm(3,3),genmmp(3,3),genswm(3,3),genswmmm(3,3),genswmmp(3,3) 572 integer :: genswp(3,3) 573 574 !************************************************************************* 575 576 !DEBUG 577 !write(std_out,*) 'symsghexa',spgroup,shubnikov,spgroupma 578 !ENDDEBUG 579 580 !The identity operation belongs to all space groups 581 symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 582 583 !Predefine some generators 584 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/) ) 585 genswmmm(:,:)=0 ; genswmmm(2,1)=-1 ; genswmmm(1,2)=-1 ; genswmmm(3,3)=-1 586 !reshape((/0,-1,0,-1,0,0,0,0,-1/),(/3,3/),(/0,0/),(/2,1/) ) 587 genswmmp(:,:)=0 ; genswmmp(2,1)=-1 ; genswmmp(1,2)=-1 ; genswmmp(3,3)=1 588 !reshape((/0,-1,0,-1,0,0,0,0,1/),(/3,3/),(/0,0/),(/2,1/) ) 589 genswp(:,:)=0 ; genswp(2,1)=1 ; genswp(1,2)=1 ; genswp(3,3)=1 590 !reshape((/0,1,0,1,0,0,0,0,1/),(/3,3/),(/0,0/),(/2,1/) ) 591 genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)=1 592 !reshape((/-1,0,0,0,-1,0,0,0,1/),(/3,3/),(/0,0/),(/2,1/) ) 593 594 !Initialize the associated translations matrix to 0 595 do ii=1,nsym 596 tnons(:,ii)= 0.0d0 597 end do 598 nogen=0 599 600 !Default non-magnetic behaviour 601 symafm(1:nsym)=1 602 603 !************************************************************************* 604 605 !Treat TRIGONAL case 606 if(143<=spgroup .and. spgroup<=167)then 607 608 ! The hexagonal axis choice (orientation) is first treated 609 if (spgaxor == 1) then 610 611 ! This matrix is common to ALL trigonal spatial groups in this orientation 612 ! (Note : this is the 3- symmetry operation) 613 symrel(:,:,2)=0 ; symrel(1,1,2)=-1 ; symrel(1,2,2)=1 ; symrel(2,1,2)=-1 ; symrel(3,3,2)=1 614 ! reshape((/-1,1,0,-1,0,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) ) 615 616 ! Assigns the generators to each space group 617 select case (spgroup) 618 case (143,146,147,148) !P3, R3, PB3, RB3 619 symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(2,2,3)=-1 ; symrel(3,3,3)=1 620 ! reshape((/0,-1,0,1,-1,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) ) 621 nogen=0 ! All symmetries have been generated 622 case (144) !P31 623 tnons(:,2)=(/0.d0,0.d0,twothird/) 624 symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(2,2,3)=-1 ; symrel(3,3,3)=1 625 ! reshape((/0,-1,0,1,-1,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) ) 626 tnons(:,3)=(/0.d0,0.d0,third/) 627 nogen=0 ! All symmetries have been generated 628 case (145) !P32 629 tnons(:,2)=(/0.d0,0.d0,third/) 630 symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(2,2,3)=-1 ; symrel(3,3,3)=1 631 ! reshape((/0,-1,0,1,-1,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) ) 632 tnons(:,3)=(/0.d0,0.d0,twothird/) 633 nogen=0 ! All symmetries have been generated 634 case (149) !P312 635 symrel(:,:,3) = genswmmm(:,:) 636 nogen=3 637 case (150,155) !P321, R32 638 symrel(:,:,3) = genswm(:,:) 639 nogen=3 640 case (151) !P3112 641 tnons(:,2)=(/0.d0,0.d0,twothird/) 642 symrel(:,:,3) = genswmmm(:,:) 643 nogen=3 644 case (152) !P3121 645 tnons(:,2)=(/0.d0,0.d0,twothird/) 646 symrel(:,:,3) = genswm(:,:) 647 nogen=3 648 case (153) !P3212 649 tnons(:,2)=(/0.d0,0.d0,third/) 650 symrel(:,:,3) = genswmmm(:,:) 651 nogen=3 652 case (154) !P3221 653 tnons(:,2)=(/0.d0,0.d0,third/) 654 symrel(:,:,3) = genswm(:,:) 655 nogen=3 656 case (156,160,164,166) !P3m1, R3m, PB3m1, RB3m 657 symrel(:,:,3) = genswmmp(:,:) 658 nogen=3 659 case (157,162) !P31m, PB31m 660 symrel(:,:,3) = genswp(:,:) 661 nogen=3 662 case (158,161,165,167) !P3c1, R3c, PB3c1, RB3c 663 symrel(:,:,3) = genswmmp(:,:) 664 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 665 nogen=3 666 case (159,163) !P31c, PB31c 667 symrel(:,:,3) = genswp(:,:) 668 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 669 nogen=3 670 end select 671 672 select case (spgroup) 673 case (146,148,155,160,166,167) 674 brvltt=7 675 end select 676 677 ! Quite simple, because the generator of even order is always the third one. 678 if(shubnikov==3)then 679 select case(spgroupma) 680 case (23,27,31,35,39,43,47,51,55,59,63,67,71,76,77,82,83,88,89,94,95,& 681 & 100,101,106,107) 682 symafm(3)=-1 683 end select 684 end if 685 686 else if (spgaxor == 2) then 687 ! The rhombohedral axis choice (orientation) is now treated 688 write(std_out,*)'rhombohedral axes' 689 ! Assignment of common three-fold rotation 690 symrel(:,:,2)=0 ; symrel(1,3,2)=1 ; symrel(3,2,2)=1 ; symrel(2,1,2)=1 691 ! reshape((/0,0,1,1,0,0,0,1,0/),(/3,3/),(/0,0/),(/2,1/) ) 692 symrel(:,:,3)=0 ; symrel(3,1,3)=1 ; symrel(2,3,3)=1 ; symrel(1,2,3)=1 693 ! reshape((/0,1,0,0,0,1,1,0,0/), (/3,3/), (/0,0/), (/2,1/) ) 694 695 select case (spgroup) 696 case (146,148) !R3 697 case (155,166) !R32, RB3m 698 symrel(:,:,4) = genswmmm(:,:) 699 nogen=4 700 case (160) !R3m 701 symrel(:,:,4) = genswp(:,:) 702 nogen=4 703 case (161,167) !R3c, RB3c 704 symrel(:,:,4) = genswp(:,:) 705 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 706 nogen=4 707 end select 708 709 if(shubnikov==3)then 710 select case(spgroupma) 711 case (47,67,71,99,101,106,107) 712 symafm(4)=-1 713 end select 714 end if 715 716 ! End selection of axis orientation 717 end if 718 719 ! End trigonal groups 720 end if 721 722 !************************************************************************* 723 724 !Treat HEXAGONAL case 725 if(168<=spgroup .and. spgroup<=194)then 726 727 ! This matrix (6) is common to most hexagonal spatial groups, except 174,187,188,189,190 728 symrel(:,:,2)=0 ; symrel(1,1,2)=1 ; symrel(1,2,2)=-1 ; symrel(2,1,2)=1 ; symrel(3,3,2)=1 729 ! reshape((/1,-1,0,1,0,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) ) 730 ! This one (6 bar) is present in the other cases 731 genm(:,:)=0 ; genm(1,2)=-1 ; genm(2,1)=1 ; genm(2,2)=-1 ; genm(3,3)=-1 732 ! reshape((/0,-1,0,1,-1,0,0,0,-1/), (/3,3/), (/0,0/), (/2,1/) ) 733 select case(spgroup) 734 case (168,175) !P6 735 nogen=2 736 case (169) !P61 737 tnons(:,2)=(/0.d0,0.d0,sixth/) 738 nogen=2 739 case (170) !P65 740 tnons(:,2)=(/0.d0,0.d0,fivesixth/) 741 nogen=2 742 case (171) !P62 743 tnons(:,2)=(/0.d0,0.d0,third/) 744 nogen=2 745 case (172) !P64 746 tnons(:,2)=(/0.d0,0.d0,twothird/) 747 nogen=2 748 case (173,176) !P63, P63/m 749 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 750 nogen=2 751 case (174) !PB6 752 symrel(:,:,2) = genm(:,:) 753 nogen=2 754 case (177) !P622 755 symrel(:,:,3) = genswm(:,:) 756 nogen=3 757 case (178) !P6122 758 tnons(:,2)=(/0.d0,0.d0,sixth/) 759 symrel(:,:,3) = genswm(:,:) 760 nogen=3 761 case (179) !P6522 762 tnons(:,2)=(/0.d0,0.d0,fivesixth/) 763 symrel(:,:,3) = genswm(:,:) 764 nogen=3 765 case (180) !P6222 766 tnons(:,2)=(/0.d0,0.d0,third/) 767 symrel(:,:,3) = genswm(:,:) 768 nogen=3 769 case (181) !P6422 770 tnons(:,2)=(/0.d0,0.d0,twothird/) 771 symrel(:,:,3) = genswm(:,:) 772 nogen=3 773 case (182) !P6322 774 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 775 symrel(:,:,3) = genswm(:,:) 776 nogen=3 777 case (183,191) !P6mm, P6/mmm 778 symrel(:,:,3) = genswp(:,:) 779 nogen=3 780 case (184,192) !P6cc, P6/mcc 781 symrel(:,:,3) = genswp(:,:) 782 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 783 nogen=3 784 case (185,193) !P63cm, P63/mcm 785 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 786 symrel(:,:,3) = genswp(:,:) 787 nogen=3 788 case (186,194) !P63mc, P63/mmc 789 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 790 symrel(:,:,3) = genswp(:,:) 791 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 792 nogen=3 793 case (187) !PB6m2 794 symrel(:,:,2)=0 ; symrel(1,2,2)=-1 ; symrel(2,1,2)=1 ; symrel(2,2,2)=-1 ; symrel(3,3,2)=-1 795 ! reshape((/0,-1,0,1,-1,0,0,0,-1/), (/3,3/), (/0,0/), (/2,1/) ) 796 symrel(:,:,3)=0 ; symrel(1,1,3)=-1 ; symrel(1,2,3)=1 ; symrel(2,2,3)=1 ; symrel(3,3,3)=1 797 ! reshape((/-1,1,0,0,1,0,0,0,1/), (/3,3/), (/0,0/), (/2,1/) ) 798 nogen=3 799 if (shubnikov==3) then 800 if (spgroupma==211) symafm(2:3)=-1 801 if (spgroupma==212) symafm(2)=-1 802 if (spgroupma==213) symafm(3)=-1 803 end if 804 case (188) !PB6c2 805 symrel(:,:,2) = genm(:,:) 806 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 807 symrel(:,:,3) = genswmmm(:,:) 808 nogen=3 809 case (189) !PB62m 810 symrel(:,:,2) = genm(:,:) 811 symrel(:,:,3) = genswp(:,:) 812 nogen=3 813 case (190) !PB62c 814 symrel(:,:,2) = genm(:,:) 815 symrel(:,:,3) = genswp(:,:) 816 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 817 nogen=3 818 end select 819 820 if(shubnikov==3)then 821 select case(spgroupma) 822 ! spgroup from 168 to 176 are OK, 177 to 194 are not done 823 case (111,115,119,123,127,131,135,139,141,145,147,152,158,164,170,& 824 & 176,182,187,193,199,205,217,224,230,237,239,247,249,257,259,267,269) 825 symafm(2)=-1 826 case(153,159,165,171,177,183,189,195,& 827 & 201,207,219,225,231,240,241,250,251,260,261,270,271) 828 symafm(3)=-1 829 case(151,157,163,169,175,181,188,194,200,206,218,223,229,236,238,246,248,256,258,266,268) 830 symafm(2:3)=-1 831 end select 832 end if 833 834 call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig) 835 836 ! End HEXAGONAL groups 837 end if 838 839 !*************************************************************************** 840 841 !DEBUG 842 !write(std_out,*) 'symsghexa : out with nogen = ',nogen 843 !ENDDEBUG 844 845 846 if (nogen>0) then 847 call bldgrp(msym,nogen,nsym,symafm,symrel,tnons) 848 end if 849 850 !DEBUG 851 !write(std_out,*)'symrel:' 852 !write(std_out,*) symrel(:,:,1:nsym) 853 !ENDDEBUG 854 855 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
SOURCE
886 subroutine symsgmono(brvltt,msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma,symafm,symrel,tnons) 887 888 implicit none 889 890 !Arguments ------------------------------------ 891 !scalars 892 integer,intent(in) :: msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma 893 integer,intent(inout) :: brvltt !vz_i 894 !arrays 895 integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i 896 real(dp),intent(inout) :: tnons(3,msym) !vz_i 897 898 !Local variables ----------------------------- 899 !scalars 900 integer :: sporder 901 character(len=1) :: brvsb 902 character(len=15) :: intsb,ptintsb,ptschsb,schsb 903 character(len=35) :: intsbl 904 !arrays 905 integer :: genmmm(3,3),genmmp(3,3),genmpm(3,3),genmpp(3,3),genpmm(3,3) 906 integer :: genpmp(3,3),genppm(3,3) 907 908 ! ************************************************************************* 909 !the identity operation belonging to all space groups 910 symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 911 912 !Predefine some generators 913 genmpp(:,:)=0 ; genmpp(1,1)=-1 ; genmpp(2,2)= 1 ; genmpp(3,3)= 1 914 genpmp(:,:)=0 ; genpmp(1,1)= 1 ; genpmp(2,2)=-1 ; genpmp(3,3)= 1 915 genppm(:,:)=0 ; genppm(1,1)= 1 ; genppm(2,2)= 1 ; genppm(3,3)=-1 916 genpmm(:,:)=0 ; genpmm(1,1)= 1 ; genpmm(2,2)=-1 ; genpmm(3,3)=-1 917 genmpm(:,:)=0 ; genmpm(1,1)=-1 ; genmpm(2,2)= 1 ; genmpm(3,3)=-1 918 genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)= 1 919 genmmm(:,:)=0 ; genmmm(1,1)=-1 ; genmmm(2,2)=-1 ; genmmm(3,3)=-1 920 921 !Default non-magnetic behaviour 922 symafm(1:nsym)=1 923 924 !assigns the generators to each space group 925 select case (spgroup) 926 case (3) ! P2 927 select case (spgaxor) 928 case (1) ! 3:b, P2_b = P2 929 symrel(:,:,2) = genmpm(:,:) 930 case (2) ! 3:a, P2_a = P2 931 symrel(:,:,2) = genpmm(:,:) 932 case (3) ! 3:c, P2_c = P2 933 symrel(:,:,2) = genmmp(:,:) 934 end select 935 if(shubnikov==3)symafm(2)=-1 936 case (4) ! P21 937 select case (spgaxor) 938 case (1) ! 3:b, P21_b = P21 939 symrel(:,:,2) = genmpm(:,:) 940 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 941 case (2) ! 3:a, P21_a = P21 942 symrel(:,:,2) = genpmm(:,:) 943 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 944 case (3) ! 3:c, P21_c = P21 945 symrel(:,:,2) = genmmp(:,:) 946 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 947 end select 948 if(shubnikov==3)symafm(2)=-1 949 case (5) ! C2 950 select case (spgaxor) 951 case (1) ! 5:b1, C2 = C2 952 symrel(:,:,2) = genmpm(:,:) 953 brvltt=4 954 case (2) ! 5:a1, B2_a = C2 955 symrel(:,:,2) = genpmm(:,:) 956 brvltt=6 957 case (3) ! 5:a2, C2_a = C2 958 symrel(:,:,2) = genpmm(:,:) 959 brvltt=4 960 case (4) ! 5:a3, I2_a = C2 961 symrel(:,:,2) = genpmm(:,:) 962 brvltt=2 963 case (5) ! 5:b2, A2_b = C2 964 symrel(:,:,2) = genmpm(:,:) 965 brvltt=5 966 case (6) ! 5:b3, I2_b = C2 967 symrel(:,:,2) = genmpm(:,:) 968 brvltt=2 969 case (7) ! 5:c1, A2_c = C2 970 symrel(:,:,2) = genmmp(:,:) 971 brvltt=5 972 case (8) ! 5:c2, B2_c = C2 973 symrel(:,:,2) = genmmp(:,:) 974 brvltt=6 975 case (9) ! 5:c3, I2_c = C2 976 symrel(:,:,2) = genmmp(:,:) 977 brvltt=2 978 end select 979 if(shubnikov==3)symafm(2)=-1 980 case (6) ! Pm 981 select case (spgaxor) 982 case (1) ! Pm_b = Pm 983 symrel(:,:,2) = genpmp(:,:) 984 case (2) ! Pm_a = Pm 985 symrel(:,:,2) = genmpp(:,:) 986 case (3) ! Pm_c = Pm 987 symrel(:,:,2) = genppm(:,:) 988 end select 989 if(shubnikov==3)symafm(2)=-1 990 case (7) ! Pc 991 select case (spgaxor) 992 case (1) ! 7:b1, Pc_b = Pc 993 symrel(:,:,2) = genpmp(:,:) 994 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 995 case (2) ! 7:a1, Pb_a = Pc 996 symrel(:,:,2) = genmpp(:,:) 997 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 998 case (3) ! 7:a2, Pn_a = Pc 999 symrel(:,:,2) = genmpp(:,:) 1000 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 1001 case (4) ! 7:a3, Pc_a = Pc 1002 symrel(:,:,2) = genmpp(:,:) 1003 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1004 case (5) ! 7:b2, Pn_b = Pc 1005 symrel(:,:,2) = genpmp(:,:) 1006 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 1007 case (6) ! 7:b3, Pa_b = Pc 1008 symrel(:,:,2) = genpmp(:,:) 1009 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1010 case (7) ! 7:c1, Pa_c = Pc 1011 symrel(:,:,2) = genppm(:,:) 1012 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1013 case (8) ! 7:c2, Pn_c = Pc 1014 symrel(:,:,2) = genppm(:,:) 1015 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 1016 case (9) ! 7:c3, Pb_c = Pb 1017 symrel(:,:,2) = genppm(:,:) 1018 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 1019 end select 1020 if(shubnikov==3)symafm(2)=-1 1021 case (8) ! Cm 1022 select case (spgaxor) 1023 case (1) ! 8:b1, Cm = Cm 1024 symrel(:,:,2) = genpmp(:,:) 1025 brvltt=4 1026 case (2) ! 8:a1, Bm_a = Cm 1027 symrel(:,:,2) = genmpp(:,:) 1028 brvltt=6 1029 case (3) ! 8:a2, Cm_a = Cm 1030 symrel(:,:,2) = genmpp(:,:) 1031 brvltt=4 1032 case (4) ! 8:a3, Im_a = Cm 1033 symrel(:,:,2) = genmpp(:,:) 1034 brvltt=2 1035 case (5) ! 8:b2, Am_b = Cm 1036 symrel(:,:,2) = genpmp(:,:) 1037 brvltt=5 1038 case (6) ! 8:b3, Im_b = Cm 1039 symrel(:,:,2) = genpmp(:,:) 1040 brvltt=2 1041 case (7) ! 8:c1, Am_c = Cm 1042 symrel(:,:,2) = genppm(:,:) 1043 brvltt=5 1044 case (8) ! 8:c2, Bm_c = Bm 1045 symrel(:,:,2) = genppm(:,:) 1046 brvltt=6 1047 case (9) ! 8:c3, Im_c = Cm 1048 symrel(:,:,2) = genppm(:,:) 1049 brvltt=2 1050 end select 1051 if(shubnikov==3)symafm(2)=-1 1052 case (9) ! Cc 1053 select case (spgaxor) 1054 case (1) ! 9:b1, Cc_b = Cc 1055 symrel(:,:,2) = genpmp(:,:) 1056 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1057 brvltt=4 1058 case (2) ! 9:a1, Bb_a = Cc 1059 symrel(:,:,2) = genmpp(:,:) 1060 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 1061 brvltt=6 1062 case (3) ! 9:a2, Cn_a = Cc 1063 symrel(:,:,2) = genmpp(:,:) 1064 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 1065 brvltt=4 1066 case (4) ! 9:a3, Ic_a = Cc 1067 symrel(:,:,2) = genmpp(:,:) 1068 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1069 brvltt=2 1070 case (5) ! 9:b2, An_b = Cc 1071 symrel(:,:,2) = genpmp(:,:) 1072 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 1073 brvltt=5 1074 case (6) ! 9:b3, Ia_b = Cc 1075 symrel(:,:,2) = genpmp(:,:) 1076 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1077 brvltt=2 1078 case (7) ! 9:c1, Aa_c = Cc 1079 symrel(:,:,2) = genppm(:,:) 1080 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1081 brvltt=5 1082 case (8) ! 9:c2, B(b+c)_c = Cc 1083 symrel(:,:,2) = genppm(:,:) 1084 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 1085 brvltt=6 1086 case (9) ! 9:c3, Ib_c = Cc 1087 symrel(:,:,2) = genppm(:,:) 1088 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 1089 brvltt=2 1090 end select 1091 if(shubnikov==3)symafm(2)=-1 1092 case (10) ! P2/m 1093 select case (spgaxor) 1094 case (1) ! 10:b, P2/m = P2/m 1095 symrel(:,:,2) = genmpm(:,:) 1096 symrel(:,:,3) = genmmm(:,:) 1097 symrel(:,:,4) = genpmp(:,:) 1098 case (2) ! 10:a, P2/m_a = P2/m 1099 symrel(:,:,2) = genpmm(:,:) 1100 symrel(:,:,3) = genmmm(:,:) 1101 symrel(:,:,4) = genmpp(:,:) 1102 case (3) ! 10:c, P2/m_c = P2/m 1103 symrel(:,:,2) = genmmp(:,:) 1104 symrel(:,:,3) = genmmm(:,:) 1105 symrel(:,:,4) = genppm(:,:) 1106 end select 1107 if(shubnikov==3)then 1108 symafm(2:4)=-1 ! Default 1109 if(spgroupma==44)symafm(4)=1 1110 if(spgroupma==45)symafm(2)=1 1111 if(spgroupma==46)symafm(3)=1 1112 end if 1113 case (11) ! P21/m 1114 select case (spgaxor) 1115 case (1) ! 11:b, P21/m = P21/m 1116 symrel(:,:,2) = genmpm(:,:) 1117 symrel(:,:,3) = genmmm(:,:) 1118 symrel(:,:,4) = genpmp(:,:) 1119 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 1120 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 1121 case (2) ! 11:a, P21/m_a = P21/m 1122 symrel(:,:,2) = genpmm(:,:) 1123 symrel(:,:,3) = genmmm(:,:) 1124 symrel(:,:,4) = genmpp(:,:) 1125 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1126 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 1127 case (3) ! 11:c, P21/m_c = P21/m 1128 symrel(:,:,2) = genmmp(:,:) 1129 symrel(:,:,3) = genmmm(:,:) 1130 symrel(:,:,4) = genppm(:,:) 1131 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1132 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 1133 end select 1134 if(shubnikov==3)then 1135 symafm(2:4)=-1 ! Default 1136 if(spgroupma==52)symafm(4)=1 1137 if(spgroupma==53)symafm(2)=1 1138 if(spgroupma==54)symafm(3)=1 1139 end if 1140 case (12) ! C2/m 1141 select case (spgaxor) 1142 case (1) ! 12:b1, C2/m = C2/m 1143 symrel(:,:,2) = genmpm(:,:) 1144 symrel(:,:,3) = genmmm(:,:) 1145 symrel(:,:,4) = genpmp(:,:) 1146 brvltt=4 1147 case (2) ! 12:a1, B2/m_a = C2/m 1148 symrel(:,:,2) = genpmm(:,:) 1149 symrel(:,:,3) = genmmm(:,:) 1150 symrel(:,:,4) = genmpp(:,:) 1151 brvltt=6 1152 case (3) ! 12:a2, C2/m_a = C2/m 1153 symrel(:,:,2) = genpmm(:,:) 1154 symrel(:,:,3) = genmmm(:,:) 1155 symrel(:,:,4) = genmpp(:,:) 1156 brvltt=4 1157 case (4) ! 12:a3, I2/m_a = C2/m 1158 symrel(:,:,2) = genpmm(:,:) 1159 symrel(:,:,3) = genmmm(:,:) 1160 symrel(:,:,4) = genmpp(:,:) 1161 brvltt=2 1162 case (5) ! 12:b2, A2/m_b = C2/m 1163 symrel(:,:,2) = genmpm(:,:) 1164 symrel(:,:,3) = genmmm(:,:) 1165 symrel(:,:,4) = genpmp(:,:) 1166 brvltt=5 1167 case (6) ! 12:b3, I2/m_b = C2/m 1168 symrel(:,:,2) = genmpm(:,:) 1169 symrel(:,:,3) = genmmm(:,:) 1170 symrel(:,:,4) = genpmp(:,:) 1171 brvltt=2 1172 case (7) ! 12:c1, A2/m_c = C2/m 1173 symrel(:,:,2) = genmmp(:,:) 1174 symrel(:,:,3) = genmmm(:,:) 1175 symrel(:,:,4) = genppm(:,:) 1176 brvltt=5 1177 case (8) ! 12:c2, B2/m_c = B2/m 1178 symrel(:,:,2) = genmmp(:,:) 1179 symrel(:,:,3) = genmmm(:,:) 1180 symrel(:,:,4) = genppm(:,:) 1181 brvltt=6 1182 case (9) ! 12:c3, I2/m_c = C2/m 1183 symrel(:,:,2) = genmmp(:,:) 1184 symrel(:,:,3) = genmmm(:,:) 1185 symrel(:,:,4) = genppm(:,:) 1186 brvltt=2 1187 end select 1188 if(shubnikov==3)then 1189 symafm(2:4)=-1 ! Default 1190 if(spgroupma==60)symafm(4)=1 1191 if(spgroupma==61)symafm(2)=1 1192 if(spgroupma==62)symafm(3)=1 1193 end if 1194 case (13) ! P2/c 1195 select case (spgaxor) 1196 case (1) ! 13:b1, P2/c = P2/c 1197 symrel(:,:,2) = genmpm(:,:) 1198 symrel(:,:,3) = genmmm(:,:) 1199 symrel(:,:,4) = genpmp(:,:) 1200 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1201 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 1202 case (2) ! 13:a1, P2/b_a = P2/c 1203 symrel(:,:,2) = genpmm(:,:) 1204 symrel(:,:,3) = genmmm(:,:) 1205 symrel(:,:,4) = genmpp(:,:) 1206 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 1207 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 1208 case (3) ! 13:a2, P2/n_a = P2/c 1209 symrel(:,:,2) = genpmm(:,:) 1210 symrel(:,:,3) = genmmm(:,:) 1211 symrel(:,:,4) = genmpp(:,:) 1212 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 1213 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 1214 case (4) ! 13:a3, P2/c_a = P2/c 1215 symrel(:,:,2) = genpmm(:,:) 1216 symrel(:,:,3) = genmmm(:,:) 1217 symrel(:,:,4) = genmpp(:,:) 1218 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1219 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 1220 case (5) ! 13:b2, P2/n_b = P2/c 1221 symrel(:,:,2) = genmpm(:,:) 1222 symrel(:,:,3) = genmmm(:,:) 1223 symrel(:,:,4) = genpmp(:,:) 1224 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 1225 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 1226 case (6) ! 13:b3, P2/a_b = P2/c 1227 symrel(:,:,2) = genmpm(:,:) 1228 symrel(:,:,3) = genmmm(:,:) 1229 symrel(:,:,4) = genpmp(:,:) 1230 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1231 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 1232 case (7) ! 13:c1, P2/a_c = P2/c 1233 symrel(:,:,2) = genmmp(:,:) 1234 symrel(:,:,3) = genmmm(:,:) 1235 symrel(:,:,4) = genppm(:,:) 1236 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1237 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 1238 case (8) ! 13:c2, P2/n_c = P2/c 1239 symrel(:,:,2) = genmmp(:,:) 1240 symrel(:,:,3) = genmmm(:,:) 1241 symrel(:,:,4) = genppm(:,:) 1242 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 1243 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 1244 case (9) ! 13:c3, P2/b_c = P2/b 1245 symrel(:,:,2) = genmmp(:,:) 1246 symrel(:,:,3) = genmmm(:,:) 1247 symrel(:,:,4) = genppm(:,:) 1248 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 1249 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 1250 end select 1251 if(shubnikov==3)then 1252 symafm(2:4)=-1 ! Default 1253 if(spgroupma==67)symafm(4)=1 1254 if(spgroupma==68)symafm(2)=1 1255 if(spgroupma==69)symafm(3)=1 1256 end if 1257 case (14) ! P21/c 1258 select case (spgaxor) 1259 case (1) ! 14:b1, P21/c_b = P21/c 1260 symrel(:,:,2) = genmpm(:,:) 1261 symrel(:,:,3) = genmmm(:,:) 1262 symrel(:,:,4) = genpmp(:,:) 1263 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 1264 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 1265 case (2) ! 14:a1, P21/a_b = P21/c 1266 symrel(:,:,2) = genpmm(:,:) 1267 symrel(:,:,3) = genmmm(:,:) 1268 symrel(:,:,4) = genmpp(:,:) 1269 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 1270 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 1271 case (3) ! 14:a2, P21/n_a = P21/c 1272 symrel(:,:,2) = genpmm(:,:) 1273 symrel(:,:,3) = genmmm(:,:) 1274 symrel(:,:,4) = genmpp(:,:) 1275 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 1276 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 1277 case (4) ! 14:a3, P21/c_a = P21/c 1278 symrel(:,:,2) = genpmm(:,:) 1279 symrel(:,:,3) = genmmm(:,:) 1280 symrel(:,:,4) = genmpp(:,:) 1281 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 1282 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 1283 case (5) ! 14:b2, P21/n_b = P21/c 1284 symrel(:,:,2) = genmpm(:,:) 1285 symrel(:,:,3) = genmmm(:,:) 1286 symrel(:,:,4) = genpmp(:,:) 1287 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 1288 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 1289 case (6) ! 14:b3, P21/a_b = P21/c 1290 symrel(:,:,2) = genmpm(:,:) 1291 symrel(:,:,3) = genmmm(:,:) 1292 symrel(:,:,4) = genpmp(:,:) 1293 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 1294 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 1295 case (7) ! 14:c1, P21/a_c = P21/c 1296 symrel(:,:,2) = genmmp(:,:) 1297 symrel(:,:,3) = genmmm(:,:) 1298 symrel(:,:,4) = genppm(:,:) 1299 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 1300 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 1301 case (8) ! 14:c2, P21/n_c = P21/c 1302 symrel(:,:,2) = genmmp(:,:) 1303 symrel(:,:,3) = genmmm(:,:) 1304 symrel(:,:,4) = genppm(:,:) 1305 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 1306 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 1307 case (9) ! 14/c3, P21/b_c = P21/b 1308 symrel(:,:,2) = genmmp(:,:) 1309 symrel(:,:,3) = genmmm(:,:) 1310 symrel(:,:,4) = genppm(:,:) 1311 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 1312 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 1313 end select 1314 if(shubnikov==3)then 1315 symafm(2:4)=-1 ! Default 1316 if(spgroupma==77)symafm(4)=1 1317 if(spgroupma==78)symafm(2)=1 1318 if(spgroupma==79)symafm(3)=1 1319 end if 1320 case (15) ! C2/c 1321 select case (spgaxor) 1322 case (1) ! 15:b1, C2/c_b = C2/c 1323 symrel(:,:,2) = genmpm(:,:) 1324 symrel(:,:,3) = genmmm(:,:) 1325 symrel(:,:,4) = genpmp(:,:) 1326 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1327 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 1328 brvltt = 4 1329 case (2) ! 15:a1, B2/b_a = C2/c 1330 symrel(:,:,2) = genpmm(:,:) 1331 symrel(:,:,3) = genmmm(:,:) 1332 symrel(:,:,4) = genmpp(:,:) 1333 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 1334 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 1335 brvltt = 6 1336 case (3) ! 15:a2, C2/n_a = C2/c 1337 symrel(:,:,2) = genpmm(:,:) 1338 symrel(:,:,3) = genmmm(:,:) 1339 symrel(:,:,4) = genmpp(:,:) 1340 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 1341 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 1342 brvltt = 4 1343 case (4) ! 15:a3, I2/c_a = C2/c 1344 symrel(:,:,2) = genpmm(:,:) 1345 symrel(:,:,3) = genmmm(:,:) 1346 symrel(:,:,4) = genmpp(:,:) 1347 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1348 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 1349 brvltt = 2 1350 case (5) ! 15:b2, A2/n_b = C2/c 1351 symrel(:,:,2) = genmpm(:,:) 1352 symrel(:,:,3) = genmmm(:,:) 1353 symrel(:,:,4) = genpmp(:,:) 1354 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 1355 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 1356 brvltt = 5 1357 case (6) ! 15:b3, I2/a_b = C2/c 1358 symrel(:,:,2) = genmpm(:,:) 1359 symrel(:,:,3) = genmmm(:,:) 1360 symrel(:,:,4) = genpmp(:,:) 1361 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1362 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 1363 brvltt = 2 1364 case (7) ! 15:c1, A2/a_c = C2/c 1365 symrel(:,:,2) = genmmp(:,:) 1366 symrel(:,:,3) = genmmm(:,:) 1367 symrel(:,:,4) = genppm(:,:) 1368 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1369 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 1370 brvltt = 5 1371 case (8) ! 15:c2, B21/b_c = C2/c 1372 symrel(:,:,2) = genmmp(:,:) 1373 symrel(:,:,3) = genmmm(:,:) 1374 symrel(:,:,4) = genppm(:,:) 1375 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 1376 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 1377 brvltt = 6 1378 case (9) ! 15:c3, I2/b_c = C2/c 1379 symrel(:,:,2) = genmmp(:,:) 1380 symrel(:,:,3) = genmmm(:,:) 1381 symrel(:,:,4) = genppm(:,:) 1382 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 1383 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 1384 brvltt = 2 1385 end select 1386 if(shubnikov==3)then 1387 symafm(2:4)=-1 ! Default 1388 if(spgroupma==87)symafm(4)=1 1389 if(spgroupma==88)symafm(2)=1 1390 if(spgroupma==89)symafm(3)=1 1391 end if 1392 end select 1393 1394 call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig) 1395 1396 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
SOURCE
1425 subroutine symsgortho(msym,nsym,shubnikov,spgaxor,spgorig,spgroup,& 1426 & spgroupma,symafm,symrel,tnons) 1427 1428 implicit none 1429 1430 !Arguments ------------------------------------ 1431 !scalars 1432 integer,intent(in) :: msym,nsym,shubnikov,spgaxor,spgorig,spgroup 1433 integer,intent(in) :: spgroupma 1434 !arrays 1435 integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i 1436 real(dp),intent(inout) :: tnons(3,msym) !vz_i 1437 1438 !Local variables ------------------------------ 1439 ! nogen = number of generators selected 1440 !scalars 1441 integer :: nogen,sporder 1442 character(len=1) :: brvsb 1443 character(len=15) :: intsb,ptintsb,ptschsb,schsb 1444 character(len=35) :: intsbl 1445 !arrays 1446 integer :: genmmm(3,3),genmmp(3,3),genmpm(3,3),genmpp(3,3),genpmm(3,3) 1447 integer :: genpmp(3,3),genppm(3,3) 1448 1449 ! ************************************************************************* 1450 1451 !DEBUG 1452 !write(std_out,*)'symsgortho ( orthorhombic groups) : enter with space group ',spgroup 1453 !ENDDEBUG 1454 1455 !The orientation of the space group: 1456 !first we will permute the input coordinates of the atoms, xred 1457 !then we will make the calculation in the "normal" space group 1458 !then the coordinates are reoriented to match the initial orientation 1459 !and finally the symrel is reoriented to correspond to the new orientation 1460 !further all the calculations are performed into the space group 1461 !with the user-defined orientation 1462 1463 !The identity operation belongs to all space groups 1464 symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 1465 1466 nogen=4 1467 1468 !Predefine some generators 1469 genmpp(:,:)=0 ; genmpp(1,1)=-1 ; genmpp(2,2)= 1 ; genmpp(3,3)= 1 1470 genpmp(:,:)=0 ; genpmp(1,1)= 1 ; genpmp(2,2)=-1 ; genpmp(3,3)= 1 1471 genppm(:,:)=0 ; genppm(1,1)= 1 ; genppm(2,2)= 1 ; genppm(3,3)=-1 1472 genpmm(:,:)=0 ; genpmm(1,1)= 1 ; genpmm(2,2)=-1 ; genpmm(3,3)=-1 1473 genmpm(:,:)=0 ; genmpm(1,1)=-1 ; genmpm(2,2)= 1 ; genmpm(3,3)=-1 1474 genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)= 1 1475 genmmm(:,:)=0 ; genmmm(1,1)=-1 ; genmmm(2,2)=-1 ; genmmm(3,3)=-1 1476 1477 1478 !For all the groups in this routine symrel(:,:,2) is the same 1479 symrel(:,:,2)=genmmp(:,:) 1480 1481 !Default non-magnetic behaviour 1482 symafm(1:nsym)=1 1483 1484 !DEBUG 1485 !write(std_out,*) 'symsgortho:',spgroup,shubnikov,spgroupma 1486 !ENDDEBUG 1487 1488 !assigns the generators to each space group 1489 select case (spgroup) 1490 ! ORTHORHOMBIC space groups 1491 case (16,21,22,23) !P222, C222, F222, I222 1492 symrel(:,:,3) = genmpm(:,:) 1493 symrel(:,:,4) = genpmm(:,:) 1494 if(shubnikov==3)symafm(3:4)=-1 1495 case (17,20) !P2221, C2221 1496 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1497 symrel(:,:,3) = genpmm(:,:) 1498 symrel(:,:,4) = genmpm(:,:) 1499 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 1500 if(shubnikov==3)then 1501 symafm(4)=-1 1502 if(spgroupma==9) symafm(3)=-1 1503 if(spgroupma==10)symafm(2)=-1 1504 if(spgroupma==33)symafm(3:4)=-1 1505 if(spgroupma==34)symafm(2)=-1 1506 end if 1507 case (18) !P21212 1508 symrel(:,:,3) = genmpm(:,:) 1509 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 1510 symrel(:,:,4) = genpmm(:,:) 1511 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 1512 if(shubnikov==3)then 1513 symafm(3)=-1 1514 if(spgroupma==18)symafm(4)=-1 1515 if(spgroupma==19)symafm(2)=-1 1516 end if 1517 case (19,24) !P212121, I212121 1518 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 1519 symrel(:,:,3) = genmpm(:,:) 1520 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 1521 symrel(:,:,4) = genpmm(:,:) 1522 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 1523 if(shubnikov==3)symafm(3:4)=-1 1524 case (25,35,38,42,44) !Pmm2, Cmm2, Amm2, Fmm2, Imm2 1525 symrel(:,:,3) = genpmp(:,:) 1526 symrel(:,:,4) = genmpp(:,:) 1527 case (26,36) !Pmc21, Cmc21 1528 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1529 symrel(:,:,3) = genpmp(:,:) 1530 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 1531 symrel(:,:,4) = genmpp(:,:) 1532 case (27,37) !Pcc2, Ccc2 1533 symrel(:,:,3) = genpmp(:,:) 1534 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 1535 symrel(:,:,4) = genmpp(:,:) 1536 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 1537 case (28,40,46) !Pma2, Ama2, Ima2 1538 symrel(:,:,3) = genpmp(:,:) 1539 tnons(:,3)=(/0.5d0,0.d0,0.d0/) 1540 symrel(:,:,4) = genmpp(:,:) 1541 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 1542 case (29) !Pca21 1543 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1544 symrel(:,:,3) = genpmp(:,:) 1545 tnons(:,3)=(/0.5d0,0.d0,0.d0/) 1546 symrel(:,:,4) = genmpp(:,:) 1547 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 1548 case (30) !Pnc2 1549 symrel(:,:,3) = genpmp(:,:) 1550 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 1551 symrel(:,:,4) = genmpp(:,:) 1552 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 1553 case (31) !Pmn21 1554 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 1555 symrel(:,:,3) = genpmp(:,:) 1556 tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 1557 symrel(:,:,4) = genmpp(:,:) 1558 case (32,41,45) !Pba2, Aba2, Iba2 1559 symrel(:,:,3) = genpmp(:,:) 1560 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 1561 symrel(:,:,4) = genmpp(:,:) 1562 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 1563 case (33) !Pna21 1564 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1565 symrel(:,:,3) = genpmp(:,:) 1566 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 1567 symrel(:,:,4) = genmpp(:,:) 1568 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 1569 case (34) !Pnn2 1570 symrel(:,:,3) = genpmp(:,:) 1571 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 1572 symrel(:,:,4) = genmpp(:,:) 1573 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 1574 case (39) !Abm2 1575 symrel(:,:,3) = genpmp(:,:) 1576 tnons(:,3)=(/0.d0,0.5d0,0.d0/) 1577 symrel(:,:,4) = genmpp(:,:) 1578 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 1579 case (43) !Fdd2 1580 symrel(:,:,3) = genpmp(:,:) 1581 tnons(:,3)=(/0.25d0,0.25d0,0.25d0/) 1582 symrel(:,:,4) = genmpp(:,:) 1583 tnons(:,4)=(/0.25d0,0.25d0,0.25d0/) 1584 case (47,65,69,71) !Pmmm, Cmmm, Fmmm, Immm 1585 symrel(:,:,3) = genmpm(:,:) 1586 symrel(:,:,4) = genpmm(:,:) 1587 case (48) !Pnnn 1588 symrel(:,:,3) = genmpm(:,:) 1589 symrel(:,:,4) = genpmm(:,:) 1590 symrel(:,:,5) = genmmm(:,:) 1591 symrel(:,:,6) = genppm(:,:) 1592 symrel(:,:,7) = genpmp(:,:) 1593 symrel(:,:,8) = genmpp(:,:) 1594 if (spgorig==1) then 1595 tnons(:,5)=(/0.5d0,0.5d0,0.5d0/) 1596 tnons(:,6)=(/0.5d0,0.5d0,0.5d0/) 1597 tnons(:,7)=(/0.5d0,0.5d0,0.5d0/) 1598 tnons(:,8)=(/0.5d0,0.5d0,0.5d0/) 1599 else 1600 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 1601 tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 1602 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 1603 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 1604 tnons(:,7)=(/0.5d0,0.d0,0.5d0/) 1605 tnons(:,8)=(/0.d0,0.5d0,0.5d0/) 1606 end if 1607 if(shubnikov==3)then 1608 if(spgroupma==259)then 1609 symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1 1610 else if(spgroupma==260)then 1611 symafm(3:4)=-1 ; symafm(7:8)=-1 1612 else if(spgroupma==261)then 1613 symafm(5:8)=-1 1614 end if 1615 end if 1616 nogen=0 1617 case (49,66) !Pccm, Cccm 1618 symrel(:,:,4) = genpmm(:,:) 1619 tnons(:,4)=(/0.d0,0.d0,0.5d0/) 1620 symrel(:,:,3) = genmpm(:,:) 1621 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 1622 case (50) !Pban 1623 symrel(:,:,3) = genmpm(:,:) 1624 symrel(:,:,4) = genpmm(:,:) 1625 symrel(:,:,5) = genmmm(:,:) 1626 symrel(:,:,6) = genppm(:,:) 1627 symrel(:,:,7) = genpmp(:,:) 1628 symrel(:,:,8) = genmpp(:,:) 1629 if (spgorig==1) then 1630 tnons(:,5)=(/0.5d0,0.5d0,0.d0/) 1631 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 1632 tnons(:,7)=(/0.5d0,0.5d0,0.d0/) 1633 tnons(:,8)=(/0.5d0,0.5d0,0.d0/) 1634 else 1635 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 1636 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 1637 tnons(:,3)=(/0.5d0,0.d0,0.d0/) 1638 tnons(:,7)=(/0.5d0,0.d0,0.d0/) 1639 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 1640 tnons(:,8)=(/0.d0,0.5d0,0.d0/) 1641 end if 1642 if(shubnikov==3)then 1643 if(spgroupma==279)then 1644 symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1 1645 else if(spgroupma==280)then 1646 symafm(3:4)=-1 ; symafm(5:6)=-1 1647 else if(spgroupma==281)then 1648 symafm(3:4)=-1 ; symafm(7:8)=-1 1649 else if(spgroupma==282)then 1650 symafm(2:8:2)=-1 1651 else if(spgroupma==283)then 1652 symafm(5:8)=-1 1653 end if 1654 end if 1655 nogen=0 1656 case (51) !Pmma 1657 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1658 symrel(:,:,3) = genmpm(:,:) 1659 symrel(:,:,4) = genpmm(:,:) 1660 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 1661 case (52) !Pnna 1662 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1663 symrel(:,:,3) = genmpm(:,:) 1664 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 1665 symrel(:,:,4) = genpmm(:,:) 1666 tnons(:,4)=(/0.d0,0.5d0,0.5d0/) 1667 case (53) !Pmna 1668 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 1669 symrel(:,:,3) = genmpm(:,:) 1670 tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 1671 symrel(:,:,4) = genpmm(:,:) 1672 case (54) !Pcca 1673 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1674 symrel(:,:,3) = genmpm(:,:) 1675 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 1676 symrel(:,:,4) = genpmm(:,:) 1677 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 1678 case (55,72) !Pbam, Ibam 1679 symrel(:,:,3) = genmpm(:,:) 1680 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 1681 symrel(:,:,4) = genpmm(:,:) 1682 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 1683 case (56) !Pccn 1684 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 1685 symrel(:,:,3) = genmpm(:,:) 1686 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 1687 symrel(:,:,4) = genpmm(:,:) 1688 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 1689 case (57) !Pbcm 1690 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1691 symrel(:,:,3) = genmpm(:,:) 1692 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 1693 symrel(:,:,4) = genpmm(:,:) 1694 tnons(:,4)=(/0.d0,0.5d0,0.d0/) 1695 case (58) !Pnnm 1696 symrel(:,:,3) = genmpm(:,:) 1697 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 1698 symrel(:,:,4) = genpmm(:,:) 1699 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 1700 case (59) !Pmmn 1701 symrel(:,:,3) = genmpm(:,:) 1702 symrel(:,:,4) = genpmm(:,:) 1703 symrel(:,:,5) = genmmm(:,:) 1704 symrel(:,:,6) = genppm(:,:) 1705 symrel(:,:,7) = genpmp(:,:) 1706 symrel(:,:,8) = genmpp(:,:) 1707 if (spgorig==1) then 1708 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 1709 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 1710 tnons(:,5)=(/0.5d0,0.5d0,0.d0/) 1711 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 1712 else 1713 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 1714 tnons(:,3)=(/0.d0,0.5d0,0.d0/) 1715 tnons(:,4)=(/0.5d0,0.d0,0.d0/) 1716 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 1717 tnons(:,7)=(/0.d0,0.5d0,0.d0/) 1718 tnons(:,8)=(/0.5d0,0.d0,0.d0/) 1719 end if 1720 if(shubnikov==3)then 1721 if(spgroupma==407)then 1722 symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1 1723 else if(spgroupma==408)then 1724 symafm(3:4)=-1 ; symafm(5:6)=-1 1725 else if(spgroupma==409)then 1726 symafm(3:4)=-1 ; symafm(7:8)=-1 1727 else if(spgroupma==410)then 1728 symafm(2:8:2)=-1 1729 else if(spgroupma==411)then 1730 symafm(5:8)=-1 1731 end if 1732 end if 1733 nogen=0 1734 case (60) !Pbcn 1735 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 1736 symrel(:,:,3) = genmpm(:,:) 1737 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 1738 symrel(:,:,4) = genpmm(:,:) 1739 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 1740 case (61,73) !Pbca, Ibca 1741 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 1742 symrel(:,:,3) = genmpm(:,:) 1743 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 1744 symrel(:,:,4) = genpmm(:,:) 1745 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 1746 case (62) !Pnma 1747 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) 1748 symrel(:,:,3) = genmpm(:,:) 1749 tnons(:,3)=(/0.d0,0.5d0,0.d0/) 1750 symrel(:,:,4) = genpmm(:,:) 1751 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 1752 case (63) !Cmcm 1753 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1754 symrel(:,:,3) = genmpm(:,:) 1755 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 1756 symrel(:,:,4) = genpmm(:,:) 1757 if(shubnikov==3)then 1758 if(spgroupma==459 .or. spgroupma==463)symafm(2:3)=-1 1759 if(spgroupma==460 .or. spgroupma==464)symafm(4)=-1 1760 if(spgroupma==460 .or. spgroupma==464)symafm(2)=-1 1761 if(spgroupma==461 .or. spgroupma==462)symafm(3:4)=-1 1762 end if 1763 case (64) !Cmca 1764 tnons(:,2)=(/0.d0,0.5d0,0.5d0/) 1765 symrel(:,:,3) = genmpm(:,:) 1766 tnons(:,3)=(/0.d0,0.5d0,0.5d0/) 1767 symrel(:,:,4) = genpmm(:,:) 1768 case (67) !Cmma 1769 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 1770 symrel(:,:,4) = genpmm(:,:) 1771 symrel(:,:,3) = genmpm(:,:) 1772 tnons(:,3)=(/0.d0,0.5d0,0.d0/) 1773 case (68) !Ccca 1774 symrel(:,:,3) = genmpm(:,:) 1775 symrel(:,:,4) = genpmm(:,:) 1776 symrel(:,:,5) = genmmm(:,:) 1777 symrel(:,:,6) = genppm(:,:) 1778 symrel(:,:,7) = genpmp(:,:) 1779 symrel(:,:,8) = genmpp(:,:) 1780 if (spgorig==1) then 1781 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) 1782 tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 1783 tnons(:,5)=(/0.d0,0.5d0,0.5d0/) 1784 tnons(:,6)=(/0.5d0,0.d0,0.5d0/) 1785 tnons(:,7)=(/0.d0,0.5d0,0.5d0/) 1786 tnons(:,8)=(/0.5d0,0.d0,0.5d0/) 1787 else 1788 tnons(:,2)=(/0.5d0,0.d0,0.d0/) 1789 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 1790 tnons(:,4)=(/0.5d0,0.d0,0.5d0/) 1791 tnons(:,6)=(/0.5d0,0.d0,0.d0/) 1792 tnons(:,7)=(/0.d0,0.d0,0.5d0/) 1793 tnons(:,8)=(/0.5d0,0.d0,0.5d0/) 1794 end if 1795 if(shubnikov==3)then 1796 if(spgroupma==513)then 1797 symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1 1798 else if(spgroupma==514)then 1799 symafm(3:4)=-1 ; symafm(5:6)=-1 1800 else if(spgroupma==515)then 1801 symafm(3:4)=-1 ; symafm(7:8)=-1 1802 else if(spgroupma==516)then 1803 symafm(2:8:2)=-1 1804 else if(spgroupma==517)then 1805 symafm(5:8)=-1 1806 end if 1807 end if 1808 nogen=0 1809 case (70) !Fddd 1810 symrel(:,:,3) = genmpm(:,:) 1811 symrel(:,:,4) = genpmm(:,:) 1812 symrel(:,:,5) = genmmm(:,:) 1813 symrel(:,:,6) = genppm(:,:) 1814 symrel(:,:,7) = genpmp(:,:) 1815 symrel(:,:,8) = genmpp(:,:) 1816 if (spgorig==1) then 1817 tnons(:,5)=(/0.25d0,0.25d0,0.25d0/) 1818 tnons(:,6)=(/0.25d0,0.25d0,0.25d0/) 1819 tnons(:,7)=(/0.25d0,0.25d0,0.25d0/) 1820 tnons(:,8)=(/0.25d0,0.25d0,0.25d0/) 1821 else 1822 tnons(:,2)=(/0.75d0,0.75d0,0.d0/) 1823 tnons(:,3)=(/0.75d0,0.d0,0.75d0/) 1824 tnons(:,4)=(/0.d0,0.75d0,0.75d0/) 1825 ! JWZ DEBUGGING BEGIN 1826 ! the following lines were present in 5.7 as of Dec 10 2008 but 1827 ! gave wrong results in a spgroup 70 spgorig 2 case (Na2SO4) 1828 ! tnons(:,5)=(/0.25d0,0.25d0,0.d0/) ! original code 1829 ! tnons(:,6)=(/0.25d0,0.d0,0.25d0/) ! original code 1830 ! tnons(:,7)=(/0.d0,0.25d0,0.25d0/) ! original code 1831 ! here are the corrected values of tnons for this case 1832 tnons(:,5)=(/0.0d0,0.0d0,0.0d0/) 1833 tnons(:,6)=(/0.25d0,0.25d0,0.0d0/) 1834 tnons(:,7)=(/0.25d0,0.0d0,0.25d0/) 1835 tnons(:,8)=(/0.0d0,0.25d0,0.25d0/) 1836 ! JWZ DEBUGGING END 1837 end if 1838 if(shubnikov==3)then 1839 if(spgroupma==529)then 1840 symafm(2:3)=-1 ; symafm(5)=-1 ; symafm(8)=-1 1841 else if(spgroupma==530)then 1842 symafm(3:4)=-1 ; symafm(7:8)=-1 1843 else if(spgroupma==531)then 1844 symafm(5:8)=-1 1845 end if 1846 end if 1847 nogen=0 1848 case (74) !Imma 1849 tnons(:,2)=(/0.d0,0.5d0,0.d0/) 1850 symrel(:,:,3) = genmpm(:,:) 1851 tnons(:,3)=(/0.d0,0.5d0,0.d0/) 1852 symrel(:,:,4) = genpmm(:,:) 1853 end select 1854 1855 if (shubnikov==3) then 1856 select case (spgroupma) 1857 case (59,68,80,89,101,113,125,137,146,158,167,174,182,189,197,& 1858 & 205,213,221,226,231,237,243,270,292,296,308,312,324,328,340,344,& 1859 & 358,370,380,384,420,424,444,448,460,464,472,476) 1860 symafm(2)=-1 1861 symafm(4)=-1 1862 case (69,90,102,114,126,147,175,190,198,206,214,244) 1863 symafm(2:3)=-1 1864 case (60,70,81,91,103,115,127,138,148,159,168,176,183,191,199,& 1865 & 207,215,222,227,232,238,245,252,268,269,293,294,309,310,& 1866 & 325,326,341,342,356,357,368,369,381,382,396,397,421,422,436,445,& 1867 & 446,461,462,473,474,484,485,494,495,504,505,524,536,& 1868 & 542,543,551,557,558) 1869 symafm(3:4)=-1 1870 case (251,267,291,295,307,311,323,327,339,343,355,367,379,383,& 1871 & 395,398,419,423,435,443,447,459,463,& 1872 & 471,475,483,486,493,496,503,506,523,535,541,544,550,556,559) 1873 symafm(2:3)=-1 1874 1875 end select 1876 end if 1877 1878 if (nogen>1)then 1879 call bldgrp(msym,nogen,nsym,symafm,symrel,tnons) 1880 end if 1881 1882 call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig) 1883 1884 !DEBUG 1885 !write(std_out,*)'symsgortho : end of symmetry assignement' 1886 !ENDDEBUG 1887 1888 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
SOURCE
1914 subroutine symsgtetra(msym,nsym,shubnikov,spgaxor,spgorig,spgroup,spgroupma,symafm,symrel,tnons) 1915 1916 implicit none 1917 1918 !Arguments ------------------------------------ 1919 !scalars 1920 integer,intent(in) :: msym,nsym,shubnikov,spgaxor,spgorig,spgroup 1921 integer,intent(in) :: spgroupma 1922 !arrays 1923 integer,intent(inout) :: symafm(msym),symrel(3,3,msym) !vz_i 1924 real(dp),intent(inout) :: tnons(3,msym) !vz_i 1925 1926 !Local variables ------------------------------ 1927 !integer :: isym 1928 !scalars 1929 integer :: nogen,sporder 1930 character(len=1) :: brvsb 1931 character(len=15) :: intsb,ptintsb,ptschsb,schsb 1932 character(len=35) :: intsbl 1933 !arrays 1934 integer :: gen4m(3,3),genmmm(3,3),genmmp(3,3),genmpm(3,3),genmpp(3,3) 1935 integer :: genpmm(3,3),genpmp(3,3),genppm(3,3) 1936 1937 ! ************************************************************************* 1938 !DEBUG 1939 !write(std_out,*) ' symsgtetra: ',spgroup,shubnikov,spgroupma 1940 !ENDDEBUG 1941 nogen=0 1942 1943 tnons(:,:)=zero 1944 !The identity operation belongs to all space groups 1945 symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 1946 1947 !The next operation belongs to most TETRAGONAL space groups, except 81,82,111,121. 1948 symrel(:,:,2)=0 ; symrel(1,2,2)=-1 ; symrel(2,1,2)=1 ; symrel(3,3,2)=1 1949 1950 !Predefine some generators 1951 genmpp(:,:)=0 ; genmpp(1,1)=-1 ; genmpp(2,2)= 1 ; genmpp(3,3)= 1 1952 genpmp(:,:)=0 ; genpmp(1,1)= 1 ; genpmp(2,2)=-1 ; genpmp(3,3)= 1 1953 genppm(:,:)=0 ; genppm(1,1)= 1 ; genppm(2,2)= 1 ; genppm(3,3)=-1 1954 genpmm(:,:)=0 ; genpmm(1,1)= 1 ; genpmm(2,2)=-1 ; genpmm(3,3)=-1 1955 genmpm(:,:)=0 ; genmpm(1,1)=-1 ; genmpm(2,2)= 1 ; genmpm(3,3)=-1 1956 genmmp(:,:)=0 ; genmmp(1,1)=-1 ; genmmp(2,2)=-1 ; genmmp(3,3)= 1 1957 genmmm(:,:)=0 ; genmmm(1,1)=-1 ; genmmm(2,2)=-1 ; genmmm(3,3)=-1 1958 gen4m(:,:)=0 ; gen4m(1,2)=-1 ; gen4m(2,1)=1 ; gen4m(3,3)=-1 1959 1960 !Default non-magnetic behaviour 1961 symafm(1:nsym)=1 1962 1963 1964 !assigns the generators to each space group 1965 select case (spgroup) 1966 ! TETRAGONAL space groups 1967 case (75,79,83,87) !P4, I4, P4/m, I4/m 1968 nogen=2 1969 case (76,80) !P41, I41 1970 tnons(:,2)=(/0.d0,0.d0,0.25d0/) 1971 nogen=2 1972 case (77,84) !P42, P42/m 1973 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 1974 nogen=2 1975 case (78) !P43 1976 tnons(:,2)=(/0.d0,0.d0,0.75d0/) 1977 nogen=2 1978 case (81,82) !PB4, IB4 1979 symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 1980 symrel(:,:,2)=0 ; symrel(1,1,2)=-1; symrel(2,2,2)=-1; symrel(3,3,2)=1 1981 symrel(:,:,3)=0 ; symrel(1,2,3)=1 ; symrel(2,1,3)=-1 ; symrel(3,3,3)=-1 1982 symrel(:,:,4)=0 ; symrel(1,2,4)=-1 ; symrel(2,1,4)=1 ; symrel(3,3,4)=-1 1983 nogen=0 1984 if (shubnikov==3) then 1985 symafm(3:4)=-1 1986 end if 1987 case (85,86,88) !P4/n 1988 symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 ; symrel(:,:,5) = -1*symrel(:,:,1) 1989 symrel(:,:,2)=0 ; symrel(1,1,2)=-1; symrel(2,2,2)=-1; symrel(3,3,2)=1 ; symrel(:,:,6) = -1*symrel(:,:,2) 1990 symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(3,3,3)=1 ; symrel(:,:,7) = -1*symrel(:,:,3) 1991 symrel(:,:,4)=0 ; symrel(1,2,4)=1 ; symrel(2,1,4)=-1 ; symrel(3,3,4)=1 ; symrel(:,:,8) = -1*symrel(:,:,4) 1992 nogen=0 1993 if(spgorig==1) then 1994 select case (spgroup) 1995 case (85) !P4/n 1996 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 1997 tnons(:,5)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 1998 case (86) !P42/n 1999 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 2000 tnons(:,5)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.5d0/) 2001 case (88) !I41/a 2002 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 2003 tnons(:,3)=(/0.0d0,0.5d0,0.25d0/) 2004 tnons(:,4)=(/0.5d0,0.0d0,0.75d0/) 2005 tnons(:,5)=(/0.0d0,0.5d0,0.25d0/) 2006 tnons(:,6)=(/0.5d0,0.0d0,0.75d0/) 2007 tnons(:,7)=(/0.5d0,0.5d0,0.5d0/) 2008 end select 2009 else 2010 select case (spgroup) 2011 case (85) !P4/n 2012 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 2013 tnons(:,3)=(/0.5d0,0.0d0,0.d0/) ; tnons(:,7)=(/0.5d0,0.0d0,0.d0/) 2014 tnons(:,4)=(/0.0d0,0.5d0,0.d0/) ; tnons(:,8)=(/0.0d0,0.5d0,0.d0/) 2015 case (86) !P42/n 2016 tnons(:,2)=(/0.5d0,0.5d0,0.0d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.0d0/) 2017 tnons(:,3)=(/0.0d0,0.5d0,0.5d0/) ; tnons(:,7)=(/0.0d0,0.5d0,0.5d0/) 2018 tnons(:,4)=(/0.5d0,0.0d0,0.5d0/) ; tnons(:,8)=(/0.5d0,0.0d0,0.5d0/) 2019 case (88) !I41/a 2020 tnons(:,2)=(/0.5d0,0.0d0,0.5d0/) ; tnons(:,6)=(/0.5d0,0.0d0,0.5d0/) 2021 tnons(:,3)=(/0.75d0,0.25d0,0.25d0/) 2022 tnons(:,7)=(/0.25d0,0.75d0,0.75d0/) 2023 tnons(:,4)=(/0.75d0,0.75d0,0.75d0/) 2024 tnons(:,8)=(/0.25d0,0.25d0,0.25d0/) 2025 end select 2026 end if 2027 if (shubnikov==3) then 2028 select case (spgroupma) 2029 case(61,69,83) 2030 symafm(3)=-1;symafm(4)=-1;symafm(7)=-1; symafm(8)=-1 2031 case(62,70,84) 2032 symafm(5)=-1;symafm(6)=-1;symafm(7)=-1; symafm(8)=-1 2033 case(63,71,85) 2034 symafm(3)=-1;symafm(5)=-1;symafm(4)=-1; symafm(6)=-1 2035 end select 2036 end if 2037 case (89,97) !P422, I422 2038 symrel(:,:,3) = genmpm(:,:) 2039 nogen=3 2040 case (90) !P4212 2041 tnons(:,2)=(/0.5d0,0.5d0,0.0d0/) 2042 symrel(:,:,3) = genmpm(:,:) 2043 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 2044 nogen=3 2045 case (91) !P4122 2046 tnons(:,2)=(/0.d0,0.d0,0.25d0/) 2047 symrel(:,:,3) = genmpm(:,:) 2048 nogen=3 2049 case (92) !P41212 2050 tnons(:,2)=(/0.5d0,0.5d0,0.25d0/) 2051 symrel(:,:,3) = genmpm(:,:) 2052 tnons(:,3)=(/0.5d0,0.5d0,0.25d0/) 2053 nogen=3 2054 case (93) !P4222 2055 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 2056 symrel(:,:,3) = genmpm(:,:) 2057 nogen=3 2058 case (94) !P42212 2059 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 2060 symrel(:,:,3) = genmpm(:,:) 2061 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 2062 nogen=3 2063 case (95) !P4322 2064 tnons(:,2)=(/0.d0,0.d0,0.75d0/) 2065 symrel(:,:,3) = genmpm(:,:) 2066 nogen=3 2067 case (96) !P43212 2068 tnons(:,2)=(/0.5d0,0.5d0,0.75d0/) 2069 symrel(:,:,3) = genmpm(:,:) 2070 tnons(:,3)=(/0.5d0,0.5d0,0.75d0/) 2071 nogen=3 2072 case (98) !I4122 2073 tnons(:,2)=(/0.d0,0.5d0,0.25d0/) 2074 symrel(:,:,3) = genmpm(:,:) 2075 tnons(:,3)=(/0.5d0,0.0d0,0.75d0/) 2076 nogen=3 2077 case (99,107,123,139) !P4mm, I4mm, P4/mmm, I4/mmm 2078 symrel(:,:,3) = genpmp(:,:) 2079 nogen=3 2080 case (100) !P4bm 2081 symrel(:,:,3) = genpmp(:,:) 2082 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 2083 nogen=3 2084 case (101,132) !P42cm, P42/mcm 2085 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 2086 symrel(:,:,3) = genpmp(:,:) 2087 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 2088 nogen=3 2089 case (102) !P42nm 2090 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 2091 symrel(:,:,3) = genpmp(:,:) 2092 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 2093 nogen=3 2094 case (103,124) !P4cc, P4/mcc 2095 symrel(:,:,3) = genpmp(:,:) 2096 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 2097 nogen=3 2098 case (104) !P4nc 2099 symrel(:,:,3) = genpmp(:,:) 2100 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 2101 nogen=3 2102 case (105,131) !P42mc, P42/mmc 2103 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 2104 symrel(:,:,3) = genpmp(:,:) 2105 nogen=3 2106 case (106,135) !P42bc, P42/mbc 2107 tnons(:,2)=(/0.d0,0.d0,0.5d0/) 2108 symrel(:,:,3) = genpmp(:,:) 2109 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 2110 nogen=3 2111 case (108) !I4cm 2112 symrel(:,:,3) = genpmp(:,:) 2113 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 2114 nogen=3 2115 case (109) !I41md 2116 tnons(:,2)=(/0.d0,0.5d0,0.25d0/) 2117 symrel(:,:,3) = genpmp(:,:) 2118 symrel(:,:,4) = genmmp(:,:) 2119 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 2120 nogen=4 2121 case (110) !I41cd 2122 tnons(:,2)=(/0.d0,0.5d0,0.25d0/) 2123 symrel(:,:,3) = genpmp(:,:) 2124 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 2125 symrel(:,:,4) = genmmp(:,:) 2126 tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 2127 nogen=4 2128 case (111,121) !PB42m, IB42m 2129 symrel(:,:,2) = gen4m(:,:) 2130 symrel(:,:,3) = genmpm(:,:) 2131 nogen=3 2132 case (112) !PB42c 2133 symrel(:,:,2) = gen4m(:,:) 2134 symrel(:,:,3) = genmpm(:,:) 2135 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 2136 nogen=3 2137 case (113) !PB421m 2138 symrel(:,:,2) = gen4m(:,:) 2139 symrel(:,:,3) = genmpm(:,:) 2140 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 2141 nogen=3 2142 case (114) !PB421c 2143 symrel(:,:,2) = gen4m(:,:) 2144 symrel(:,:,3) = genmpm(:,:) 2145 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 2146 nogen=3 2147 case (115,119) !PB4m2,IB4m2 2148 symrel(:,:,2) = gen4m(:,:) 2149 symrel(:,:,3) = genpmp(:,:) 2150 nogen=3 2151 case (116,120) !PB4c2, IB4c2 2152 symrel(:,:,2) = gen4m(:,:) 2153 symrel(:,:,3) = genpmp(:,:) 2154 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 2155 nogen=3 2156 case (117) !PB4b2 2157 symrel(:,:,2) = gen4m(:,:) 2158 symrel(:,:,3) = genpmp(:,:) 2159 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 2160 nogen=3 2161 case (118) !PB4n2 2162 symrel(:,:,2) = gen4m(:,:) 2163 symrel(:,:,3) = genpmp(:,:) 2164 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 2165 nogen=3 2166 case (122) !IB42d 2167 symrel(:,:,2)=0 ; symrel(1,2,2)=-1 ; symrel(2,1,2)=1 ; symrel(3,3,2)=-1 2168 symrel(:,:,3) = genmpm(:,:) 2169 tnons(:,3)=(/0.5d0,0.d0,0.75d0/) 2170 nogen=3 2171 case (125,126,129,130,133,134,137,138,141,142) 2172 symrel(:,:,1)=0 ; symrel(1,1,1)=1 ; symrel(2,2,1)=1 ; symrel(3,3,1)=1 2173 symrel(:,:,2)=0 ; symrel(1,1,2)=-1; symrel(2,2,2)=-1; symrel(3,3,2)=1 2174 symrel(:,:,3)=0 ; symrel(1,2,3)=-1 ; symrel(2,1,3)=1 ; symrel(3,3,3)=1 2175 symrel(:,:,4)=0 ; symrel(1,2,4)=1 ; symrel(2,1,4)=-1 ; symrel(3,3,4)=1 2176 symrel(:,:,5) = genmpm(:,:) 2177 symrel(:,:,6) = genmmm(:,:) 2178 if (spgorig==1) then 2179 select case (spgroup) 2180 case (125) !P4/nbm 2181 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 2182 case (126) !P4/nnc 2183 tnons(:,6)=(/0.5d0,0.5d0,0.5d0/) 2184 case (129) !P4/nmm 2185 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 2186 tnons(:,5)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 2187 case (130) !P4/ncc 2188 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.d0/) 2189 tnons(:,5)=(/0.5d0,0.5d0,0.5d0/) 2190 tnons(:,6)=(/0.5d0,0.5d0,0.d0/) 2191 case (133) !P42/nbc 2192 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 2193 tnons(:,5)=(/0.d0,0.d0,0.5d0/) 2194 tnons(:,6)=(/0.5d0,0.5d0,0.5d0/) 2195 case (134) !P42/nnm 2196 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 2197 tnons(:,6)=(/0.5d0,0.5d0,0.5d0/) 2198 case (137) !P42/nmc 2199 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 2200 tnons(:,5)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.5d0/) 2201 case (138) !P42/ncm 2202 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) ; tnons(:,4)=(/0.5d0,0.5d0,0.5d0/) 2203 tnons(:,5)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,6)=(/0.5d0,0.5d0,0.5d0/) 2204 case (141) !I41/amd 2205 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 2206 tnons(:,3)=(/0.d0,0.5d0,0.25d0/) 2207 tnons(:,4)=(/0.5d0,0.d0,0.75d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.75d0/) 2208 tnons(:,6)=(/0.d0,0.5d0,0.25d0/) 2209 case (142) !I41/acd 2210 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 2211 tnons(:,3)=(/0.d0,0.5d0,0.25d0/) 2212 tnons(:,4)=(/0.5d0,0.d0,0.75d0/) 2213 tnons(:,5)=(/0.5d0,0.d0,0.25d0/) 2214 tnons(:,6)=(/0.d0,0.5d0,0.25d0/) 2215 end select 2216 else 2217 select case (spgroup) 2218 case (125) !P4/nbm 2219 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.d0/) 2220 tnons(:,4)=(/0.0d0,0.5d0,0.d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.d0/) 2221 case (126) !P4/nnc 2222 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.d0/) 2223 tnons(:,4)=(/0.0d0,0.5d0,0.d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.5d0/) 2224 case (129) !P4/nmm 2225 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.d0/) 2226 tnons(:,4)=(/0.0d0,0.5d0,0.d0/) ; tnons(:,5)=(/0.0d0,0.5d0,0.d0/) 2227 case (130) !P4/ncc 2228 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.d0/) 2229 tnons(:,4)=(/0.0d0,0.5d0,0.d0/) ; tnons(:,5)=(/0.0d0,0.5d0,0.5d0/) 2230 case (133) !P42/nbc 2231 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 2232 tnons(:,4)=(/0.0d0,0.5d0,0.5d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.d0/) 2233 case (134) !P42/nnm 2234 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 2235 tnons(:,4)=(/0.0d0,0.5d0,0.5d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.5d0/) 2236 case (137) !P42/nmc 2237 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 2238 tnons(:,4)=(/0.0d0,0.5d0,0.5d0/) ; tnons(:,5)=(/0.d0,0.5d0,0.d0/) 2239 case (138) !P42/ncm 2240 tnons(:,2)=(/0.5d0,0.5d0,0.d0/) ; tnons(:,3)=(/0.5d0,0.d0,0.5d0/) 2241 tnons(:,4)=(/0.0d0,0.5d0,0.5d0/) ; tnons(:,5)=(/0.d0,0.5d0,0.5d0/) 2242 case (141) !I41/amd 2243 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) ; tnons(:,3)=(/0.25d0,0.75d0,0.25d0/) 2244 tnons(:,4)=(/0.25d0,0.25d0,0.75d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.5d0/) 2245 case (142) !I41/acd 2246 tnons(:,2)=(/0.5d0,0.d0,0.5d0/) ; tnons(:,3)=(/0.25d0,0.75d0,0.25d0/) 2247 tnons(:,4)=(/0.25d0,0.25d0,0.75d0/) ; tnons(:,5)=(/0.5d0,0.d0,0.d0/) 2248 end select 2249 end if 2250 nogen=6 2251 case (127) !P4/mbm 2252 symrel(:,:,3) = genmpm(:,:) 2253 tnons(:,3)=(/0.5d0,0.5d0,0.d0/) 2254 nogen=3 2255 case (128) !P4/mnc 2256 symrel(:,:,3) = genmpm(:,:) 2257 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 2258 nogen=3 2259 case (136) !P42/mnm 2260 tnons(:,2)=(/0.5d0,0.5d0,0.5d0/) 2261 symrel(:,:,3) = genmpm(:,:) 2262 tnons(:,3)=(/0.5d0,0.5d0,0.5d0/) 2263 nogen=3 2264 case (140) !I4/mcm 2265 symrel(:,:,3) = genmpm(:,:) 2266 tnons(:,3)=(/0.d0,0.d0,0.5d0/) 2267 nogen=3 2268 end select 2269 2270 if(shubnikov==3)then 2271 select case(spgroupma) 2272 case(3,9,15,21,27,31,45,47,53,55,77,79,& 2273 & 89,97,105,113,121,129,137,145,153,159,166,174,182,190,198,206,& 2274 & 214,222,230,236,242,248,254,262,270,278,285,293,301,309,317,& 2275 & 323,330,336,343,346,355,358,391,392,403,404,& 2276 & 439,442,451,454,487,490,499,500,535,& 2277 & 538,545,546) 2278 symafm(2)=-1 2279 case(90,98,106,114,122,130,138,146,154,160,167,175,183,191,199,& 2280 & 207,215,223,231,237,243,249,255,263,271,279,287,295,303,311,319,& 2281 & 325,331,337,345,347,357,359,389,393,401,405,441,443,453,455,& 2282 & 489,491,497,501,537,539,543,547) 2283 symafm(3)=-1 2284 case(91,99,107,115,123,131,139,147,155,161,165,173,181,189,& 2285 & 197,205,213,221,229,235,241,247,253,261,269,277,286,294,302,310,& 2286 & 318,324,329,335,342,344,354,356,390,394,402,406,& 2287 & 438,440,450,452,486,488,498,502,534,536,544,548) 2288 symafm(2:3)=-1 2289 ! case(365,377,413,425,461,473,509,521) !XG230719 2290 ! symafm(2)=-1 2291 ! symafm(5:6)=-1 2292 case(366,378,414,426,462,474,510,522,554,564) 2293 symafm(3:5)=-1 2294 case(367,379,415,427,463,475,511,523,555,565) 2295 symafm(3:4)=-1 2296 case(368,380,416,428,464,476,512,524,556,566) 2297 symafm(3:4)=-1 2298 symafm(6)=-1 2299 ! case(369,381,417,429,465,477,513,525) !XG230719 2300 ! symafm(2)=-1 2301 ! symafm(5)=-1 2302 case(370,382,418,430,466,478,514,526,558,568) 2303 symafm(3:6)=-1 2304 case(371,383,419,431,467,479,515,527,559,569) 2305 symafm(6)=-1 2306 ! case(553,563) 2307 case(365,377,413,425,461,473,509,521,553,563) !XG230719 2308 symafm(5:6)=-1 2309 ! case(557,567) 2310 case(369,381,417,429,465,477,513,525,557,567) !XG230719 2311 symafm(5)=-1 2312 end select 2313 end if 2314 2315 !DEBUG 2316 !write(std_out,*)' symsgtetra : symrel(:,:,6)=',symrel(:,:,6) 2317 !ENDDEBUG 2318 2319 if (nogen>1) then 2320 call bldgrp(msym,nogen,nsym,symafm,symrel,tnons) 2321 end if 2322 2323 call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig) 2324 2325 !DEBUG 2326 !write(std_out,*)' symsgtetra : exit' 2327 !do isym=1,nsym 2328 !write(std_out,'(i3,2x,9i3,3es13.3,i3)') isym,symrel(:,:,isym),tnons(:,isym),symafm(isym) 2329 !end do 2330 !ENDDEBUG 2331 2332 end subroutine symsgtetra