Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 : ! This file is part of FLEUR and available as free software under the conditions
4 : ! of the MIT license as expressed in the LICENSE file in more detail.
5 : !--------------------------------------------------------------------------------
6 :
7 : MODULE m_grp_k
8 : USE m_juDFT
9 :
10 : IMPLICIT NONE
11 :
12 : PRIVATE
13 : PUBLIC :: grp_k, euler
14 :
15 : CONTAINS
16 :
17 0 : SUBROUTINE grp_k(sym,mrot_k,cell,bk,nclass,nirr,char_table,&
18 0 : grpname,irrname,su,sp_alph,sp_beta)
19 :
20 : !************************************************************
21 : !
22 : ! Determines the group of k, returns the number of classes, the name of the
23 : ! group, the irreducible representations and the character table.
24 : ! All the groups are not implemented yet, and the identification is not tested for
25 : ! the groups.
26 : !
27 : ! The subroutine works also with double groups, however, no double groups are
28 : ! tabulated at the moment.
29 : !
30 : !
31 : ! Character tables are taken from
32 : ! T. Inui, Y Tanabe, and Y. Onodera, "Group theory and its applications in physics",
33 : ! Springer (1996)
34 : !
35 : ! Jussi Enkovaara, Juelich 2004
36 : !**************************************************************
37 :
38 : ! USE m_mrot2su
39 : USE m_inv3
40 : USE m_constants
41 : USE m_socsym, ONLY : soc_sym, cross
42 : USE m_types
43 : IMPLICIT NONE
44 :
45 : TYPE(t_sym),INTENT(IN) :: sym
46 : TYPE(t_cell),INTENT(IN) :: cell
47 :
48 :
49 : REAL,INTENT(IN) :: bk(3)
50 : COMPLEX, INTENT(OUT) :: char_table(:,:)
51 : INTEGER, INTENT(OUT) :: mrot_k(:,:,:),nclass,nirr
52 : CHARACTER(LEN=7), INTENT(OUT) :: grpname
53 : CHARACTER(LEN=5), INTENT(OUT) :: irrname(:)
54 : COMPLEX, OPTIONAL, INTENT(OUT) :: su(:,:,:)
55 : REAL, OPTIONAL, INTENT(IN) :: sp_alph,sp_beta ! spin quant. axis
56 :
57 : ! locals
58 : INTEGER :: nopk,c
59 0 : LOGICAL, ALLOCATABLE :: error(:)
60 : INTEGER :: mrot2(3,3,48)
61 : REAL :: ktest(3),rtmp,rtmp2(4),kt(3)
62 : INTEGER :: mtest(3,3),mtmp(3,3),mtmpinv(3,3),munit(3,3)
63 : INTEGER :: n,n2 ,i,itmp
64 : INTEGER :: elem(48),belongs(48),members(48)
65 : LOGICAL :: soc,l_sorted
66 : ! COMPLEX :: sutmp(2,2), sutmpinv(2,2), su2(2,2,48)
67 : COMPLEX :: sutest(2,2)
68 : INTEGER :: d,det(48),rot(12),rot_type(48)
69 : REAL :: alpha,beta,gamma,theta(48),rax(4,48)
70 : REAL,PARAMETER :: eps=0.000001
71 : COMPLEX :: omega
72 : COMPLEX,PARAMETER :: one=CMPLX(1.0,0.0),zero=CMPLX(0.0,0.0)
73 : COMPLEX,PARAMETER :: imi=CMPLX(0.0,1.0)
74 :
75 :
76 0 : soc=.FALSE.
77 0 : IF (PRESENT(su)) soc=.TRUE.
78 : !sym%nop=SIZE(sym%mrot,3)
79 :
80 0 : ALLOCATE(error(sym%nop))
81 0 : error=.FALSE.
82 : ! Reduce the symmetry due to spin-orbit
83 0 : IF (soc.AND.PRESENT(sp_alph)) THEN
84 0 : CALL soc_sym(sym%nop,sym%mrot,sp_beta,sp_alph,cell%amat,error)!keep
85 : ENDIF
86 :
87 : ! determine the group of k
88 : nopk=0
89 0 : ksymloop: DO n=1,sym%nop
90 0 : ktest(1)=bk(1)-sym%mrot(1,1,n)*bk(1)-sym%mrot(2,1,n)*bk(2)-sym%mrot(3,1,n)*bk(3)
91 0 : ktest(2)=bk(2)-sym%mrot(1,2,n)*bk(1)-sym%mrot(2,2,n)*bk(2)-sym%mrot(3,2,n)*bk(3)
92 0 : ktest(3)=bk(3)-sym%mrot(1,3,n)*bk(1)-sym%mrot(2,3,n)*bk(2)-sym%mrot(3,3,n)*bk(3)
93 : IF (( ABS( ktest(1) - NINT(ktest(1)) ) < eps ) .AND.&
94 : & ( ABS( ktest(2) - NINT(ktest(2)) ) < eps ) .AND.&
95 0 : & ( ABS( ktest(3) - NINT(ktest(3)) ) < eps ) .AND.&
96 0 : & (.NOT.error(n)) ) THEN
97 0 : nopk=nopk+1
98 0 : mrot_k(:,:,nopk)=sym%mrot(:,:,n)
99 : CYCLE ksymloop
100 : ENDIF
101 : ENDDO ksymloop
102 :
103 0 : DEALLOCATE(error)
104 :
105 : ! Determine the spin-rotations
106 : ! Double groups not used at the moment, the groups are classified
107 : ! without the spin rotations
108 : ! IF (soc) CALL mrot2su(mrot_k(:,:,1:nopk),amat,su)
109 :
110 : ! identify the group
111 : ! first determine the classes
112 0 : members=1
113 0 : nclass=1
114 0 : elem(nclass)=1
115 0 : belongs(1)=1
116 : ! IF (soc) THEN
117 : ! double group
118 : ! DO n=1,nopk
119 : ! mrot_k(:,:,n+nopk)=mrot_k(:,:,n)
120 : ! su(:,:,n+nopk)=-su(:,:,n)
121 : ! ENDDO
122 : ! nopk=2*nopk
123 : ! ENDIF
124 :
125 0 : cloop: DO n=2,nopk
126 0 : DO n2=1,nopk
127 0 : mtmp=mrot_k(:,:,n2)
128 0 : CALL inv3(mrot_k(:,:,n2),mtmpinv,d)
129 0 : mtest=MATMUL(MATMUL(mtmp,mrot_k(:,:,n)),mtmpinv)
130 : ! IF (soc) THEN
131 : ! sutmp=su(:,:,n2)
132 : ! sutmpinv=CONJG(TRANSPOSE(sutmp))
133 : ! sutest=MATMUL(MATMUL(sutmp,su(:,:,n)),sutmpinv)
134 : ! ENDIF
135 0 : DO c=1,nclass
136 : ! IF (soc) THEN
137 : ! IF (ALL((mtest-mrot_k(:,:,elem(c))).EQ.0).AND.
138 : !c & ALL(ABS(REAL(sutest-su(:,:,elem(c)))).LE.0.0001).AND.
139 : ! & ALL(ABS(AIMAG(sutest-su(:,:,elem(c)))).LE.0.0001)) THEN
140 : ! belongs(n)=c
141 : ! CYCLE classloop
142 : ! ENDIF
143 : ! ELSE
144 0 : IF (ALL((mtest-mrot_k(:,:,elem(c))).EQ.0)) THEN
145 0 : belongs(n)=c
146 0 : members(c)=members(c)+1
147 0 : CYCLE cloop
148 : ENDIF
149 : ! ENDIF
150 : ENDDO
151 : ENDDO
152 0 : nclass=nclass+1
153 0 : elem(nclass)=n
154 0 : belongs(n)=nclass
155 : ENDDO cloop
156 0 : nirr=nclass
157 0 : IF (soc) nirr=2*nclass
158 :
159 : ! sort the classes according the rotation type
160 0 : mrot2(:,:,1:nopk)=mrot_k(:,:,1:nopk)
161 0 : DO c=1,nclass
162 0 : mrot_k(:,:,c)=mrot2(:,:,elem(c))
163 : ENDDO
164 : ! IF (soc) THEN
165 : ! su2(:,:,1:nopk)=su(:,:,1:nopk)
166 : ! DO c=1,nclass
167 : ! su(:,:,c)=su2(:,:,elem(c))
168 : ! ENDDO
169 : ! ENDIF
170 :
171 :
172 :
173 : ! identify the group
174 0 : rot=0
175 0 : rot_type=0
176 : ! rot_type: 2=two fold rot, 3=three fold rot, etc, 7=identity
177 : ! 8=two fold improper rot, ...
178 0 : munit=0
179 0 : munit(1,1)=1
180 0 : munit(2,2)=1
181 0 : munit(3,3)=1
182 : ! determine the number of different rotations
183 0 : DO c=1,nclass
184 : det(c)=&
185 : & mrot_k(1,1,c)*mrot_k(2,2,c)*mrot_k(3,3,c)+&
186 : & mrot_k(1,2,c)*mrot_k(2,3,c)*mrot_k(3,1,c)+&
187 : & mrot_k(2,1,c)*mrot_k(3,2,c)*mrot_k(1,3,c)-&
188 : & mrot_k(1,3,c)*mrot_k(2,2,c)*mrot_k(3,1,c)-&
189 : & mrot_k(2,3,c)*mrot_k(3,2,c)*mrot_k(1,1,c)-&
190 0 : & mrot_k(2,1,c)*mrot_k(1,2,c)*mrot_k(3,3,c)
191 0 : mtest=det(c)*mrot_k(:,:,c)
192 0 : rotloop: DO i=1,6
193 0 : IF (ALL((mtest-munit).EQ.0)) THEN
194 0 : rot(i+(1-det(c))*3)=rot(i+(1-det(c))*3)+1
195 0 : IF (i.EQ.1) THEN
196 0 : rot_type(c)=7-(1-det(c))*3
197 : ELSE
198 0 : rot_type(c)=i-(1-det(c))*3
199 : ENDIF
200 0 : EXIT rotloop
201 : ENDIF
202 0 : mtest=MATMUL(det(c)*mrot_k(:,:,c),mtest)
203 : ENDDO rotloop
204 0 : IF (ANY((mtest-munit).NE.0)) THEN
205 0 : WRITE(oUnit,*) 'grp_k: Cannot find the order of rotation'
206 : ENDIF
207 0 : IF ((rot(5).GT.0).OR.(rot(11).GT.0)) THEN
208 0 : WRITE(oUnit,*) 'grp_k: 5 fold rotation found!'
209 : ENDIF
210 0 : CALL euler(c,sym,cell,alpha,beta,gamma)
211 0 : CALL rotaxis(alpha,beta,gamma,rax(1:4,c),theta(c))
212 0 : IF (soc) THEN
213 0 : su(1,1,c)=COS(beta/2.0)*EXP(CMPLX(0.0,-(alpha+gamma)/2.0))
214 0 : su(1,2,c)=-SIN(beta/2.0)*EXP(CMPLX(0.0,-(alpha-gamma)/2.0))
215 0 : su(2,1,c)=SIN(beta/2.0)*EXP(CMPLX(0.0,(alpha-gamma)/2.0))
216 0 : su(2,2,c)=COS(beta/2.0)*EXP(CMPLX(0.0,(alpha+gamma)/2.0))
217 : ENDIF
218 : ENDDO
219 :
220 : !<-- Sort the classes
221 : ! The group elements are sorted in the following way:
222 : ! First the proper rotations, then improper ones, with increasing rotation angle
223 : ! Rotations with the same angle are arranged with increasing magnitude of the rotation
224 : ! axis, e.g. (I, 90 deg around 001, 90 around 110, mirror, ...
225 :
226 : l_sorted=.FALSE.
227 0 : DO WHILE (.NOT.l_sorted)
228 0 : l_sorted=.TRUE.
229 0 : DO c=1,nclass-1
230 0 : IF (rot_type(c).LT.rot_type(c+1)) THEN
231 0 : mtest=mrot_k(:,:,c)
232 0 : mrot_k(:,:,c)=mrot_k(:,:,c+1)
233 0 : mrot_k(:,:,c+1)=mtest
234 0 : rtmp=theta(c)
235 0 : theta(c)=theta(c+1)
236 0 : theta(c+1)=rtmp
237 0 : rtmp2=rax(:,c)
238 0 : rax(:,c)=rax(:,c+1)
239 0 : rax(:,c+1)=rtmp2
240 0 : itmp=det(c)
241 0 : det(c)=det(c+1)
242 0 : det(c+1)=itmp
243 0 : itmp=rot_type(c)
244 0 : rot_type(c)=rot_type(c+1)
245 0 : rot_type(c+1)=itmp
246 0 : itmp=members(c)
247 0 : members(c)=members(c+1)
248 0 : members(c+1)=itmp
249 0 : IF (soc) THEN
250 0 : sutest=su(:,:,c)
251 0 : su(:,:,c)=su(:,:,c+1)
252 0 : su(:,:,c+1)=sutest
253 : ENDIF
254 : l_sorted=.FALSE.
255 : ENDIF
256 0 : IF ((rot_type(c).EQ.rot_type(c+1)).AND.&
257 : & (theta(c).GT.theta(c+1))) THEN
258 : ! IF (theta(c)+(1-det(c))*2.0*pi.GT.
259 : ! & theta(c+1)+(1-det(c+1))*2.0*pi) THEN
260 0 : mtest=mrot_k(:,:,c)
261 0 : mrot_k(:,:,c)=mrot_k(:,:,c+1)
262 0 : mrot_k(:,:,c+1)=mtest
263 0 : rtmp=theta(c)
264 0 : theta(c)=theta(c+1)
265 0 : theta(c+1)=rtmp
266 0 : rtmp2=rax(:,c)
267 0 : rax(:,c)=rax(:,c+1)
268 0 : rax(:,c+1)=rtmp2
269 0 : itmp=det(c)
270 0 : det(c)=det(c+1)
271 0 : det(c+1)=itmp
272 0 : itmp=members(c)
273 0 : members(c)=members(c+1)
274 0 : members(c+1)=itmp
275 0 : IF (soc) THEN
276 0 : sutest=su(:,:,c)
277 0 : su(:,:,c)=su(:,:,c+1)
278 0 : su(:,:,c+1)=sutest
279 : ENDIF
280 : l_sorted=.FALSE.
281 : ENDIF
282 : IF ((rot_type(c).EQ.rot_type(c+1)).AND.&
283 0 : & (ABS(theta(c)-theta(c+1)).LT.0.0001).AND. &
284 0 : & (rax(4,c).GT.rax(4,c+1))) THEN
285 : ! IF ((ABS(theta(c)+(1-det(c))*2.0*pi-
286 : ! & theta(c+1)-(1-det(c+1))*2.0*pi).LT.0.0001).AND.
287 : ! & (rax(4,c).GT.rax(4,c+1))) THEN
288 0 : mtest=mrot_k(:,:,c)
289 0 : mrot_k(:,:,c)=mrot_k(:,:,c+1)
290 0 : mrot_k(:,:,c+1)=mtest
291 0 : rtmp2=rax(:,c)
292 0 : rax(:,c)=rax(:,c+1)
293 0 : rax(:,c+1)=rtmp2
294 0 : itmp=members(c)
295 0 : members(c)=members(c+1)
296 0 : members(c+1)=itmp
297 0 : IF (soc) THEN
298 0 : sutest=su(:,:,c)
299 0 : su(:,:,c)=su(:,:,c+1)
300 0 : su(:,:,c+1)=sutest
301 : ENDIF
302 : l_sorted=.FALSE.
303 : ENDIF
304 : ENDDO
305 : ENDDO
306 : !>
307 :
308 :
309 0 : WRITE(24,110) bk
310 0 : DO c=1,nclass
311 0 : IF (det(c).EQ.1) THEN
312 0 : WRITE(24,111) NINT(theta(c)*180/pi_const),(rax(1:3,c)),members(c)
313 : ELSE
314 0 : WRITE(24,112) NINT(theta(c)*180/pi_const),(rax(1:3,c)),members(c)
315 : ENDIF
316 : ENDDO
317 : 110 FORMAT('Symmetry operations of group of k=',3f6.3)
318 : 111 FORMAT(i3,1x,'degree proper rotation around ',3f6.2,3x,&
319 : & i3,' members in class')
320 : 112 FORMAT(i3,1x,'degree improper rotation around ',3f6.2,3x,i3&
321 : & ,' members in class')
322 :
323 :
324 : !<-- Character tables
325 :
326 0 : char_table=0.0
327 0 : char_table(1,:)=1.0
328 0 : grpname='Unknown'
329 0 : irrname='Unkno' ! Only 5 characters wide, cf. ../sympsi.F.
330 : ! First look the number of classes, within groups with the same number of classes
331 : ! check the number of different rotations
332 0 : SELECT CASE(nclass)
333 : CASE(1) ! only C1
334 0 : grpname='C1'
335 0 : char_table(1,1)=1.0
336 0 : irrname(1)='Gam1'
337 0 : IF (soc) THEN
338 0 : char_table(2,1)=1.0
339 0 : irrname(2)='Gam2'
340 : ENDIF
341 : CASE(2) ! C-1, C2 and Cs
342 0 : IF (rot(2).GT.0) THEN
343 0 : grpname='C2'
344 0 : char_table(1,1:2)=(/1.0, 1.0/)
345 0 : char_table(2,1:2)=(/1.0, -1.0/)
346 0 : irrname(1)='Gam1'
347 0 : irrname(2)='Gam2'
348 : ELSE
349 0 : grpname='C1h'
350 0 : char_table(1,1:2)=(/1.0, 1.0/)
351 0 : char_table(2,1:2)=(/1.0, -1.0/)
352 0 : irrname(1)='Gam1'
353 0 : irrname(2)='Gam2'
354 0 : IF (soc) THEN
355 0 : char_table(3,1:2)=(/one, -imi/)
356 0 : char_table(4,1:2)=(/one, imi/)
357 0 : irrname(3)='Gam3'
358 0 : irrname(4)='Gam4'
359 : ENDIF
360 : ENDIF
361 : CASE(3) ! C3, D3 and C3v
362 0 : IF (ANY(det(1:3).EQ.-1)) THEN
363 0 : grpname='C3v'
364 0 : char_table(1,1:3)=(/1.0, 1.0, 1.0/)
365 0 : char_table(2,1:3)=(/1.0, 1.0, -1.0/)
366 0 : char_table(3,1:3)=(/2.0, -1.0, 0.0/)
367 0 : irrname(1)='Lam1'
368 0 : irrname(2)='Lam2'
369 0 : irrname(3)='Lam3'
370 0 : IF (soc) THEN
371 0 : nirr=6
372 0 : char_table(4,1:3)=(/2.0, 1.0, 0.0/)
373 0 : char_table(5,1:3)=(/one, -one, imi/)
374 0 : char_table(6,1:3)=(/one, -one, -imi/)
375 0 : irrname(4)='Lam6'
376 0 : irrname(5)='Lam4'
377 0 : irrname(6)='Lam5'
378 : ENDIF
379 : ELSE
380 0 : grpname='C3'
381 0 : omega=EXP(CMPLX(0.0,-2*pi_const/3.0))
382 0 : char_table(1,1:3)=(/1.0, 1.0, 1.0/)
383 0 : char_table(2,1:3)=(/one, omega, omega**2/)
384 0 : char_table(3,1:3)=(/one, omega**2, omega /)
385 0 : irrname(1)='Gam1'
386 0 : irrname(2)='Gam2'
387 0 : irrname(3)='Gam3'
388 : ENDIF
389 : CASE(4) ! C2h, D2, C2v, C4, S4, and T
390 0 : IF((rot(2).EQ.1).AND.(rot(8).EQ.2)) THEN
391 0 : grpname='C2v'
392 0 : char_table(1,1:4)=(/1.0, 1.0, 1.0, 1.0/)
393 0 : char_table(2,1:4)=(/1.0, 1.0, -1.0, -1.0/)
394 0 : char_table(3,1:4)=(/1.0, -1.0, -1.0, 1.0/)
395 0 : char_table(4,1:4)=(/1.0, -1.0, 1.0, -1.0/)
396 0 : irrname(1)='Z1'
397 0 : irrname(2)='Z2'
398 0 : irrname(3)='Z3'
399 0 : irrname(4)='Z4'
400 0 : ELSE IF (rot(3).GT.1) THEN
401 0 : grpname='T'
402 0 : omega=EXP(CMPLX(0.0,-2*pi_const/3.0))
403 0 : char_table(1,1:4)=(/1.0, 1.0, 1.0, 1.0/)
404 0 : char_table(2,1:4)=(/one, one, omega, omega**2/)
405 0 : char_table(3,1:4)=(/one, one, omega**2, omega/)
406 0 : char_table(4,1:4)=(/3.0, -1.0, 0.0, 0.0/)
407 0 : irrname(1)='Gam1'
408 0 : irrname(2)='Gam2'
409 0 : irrname(3)='Gam3'
410 0 : irrname(4)='Gam4'
411 0 : ELSE IF ((rot(4).GT.1).OR.(rot(10).GT.1)) THEN
412 0 : IF (rot(4).GT.1) grpname='C4'
413 0 : IF (rot(10).GT.1) grpname='S4'
414 0 : char_table(1,1:4)=(/1.0, 1.0, 1.0, 1.0/)
415 0 : char_table(2,1:4)=(/1.0, -1.0, 1.0, -1.0/)
416 0 : char_table(3,1:4)=(/one, -imi, -one, imi/)
417 0 : char_table(4,1:4)=(/one, imi, -one, -imi/)
418 0 : irrname(1)='Gam1'
419 0 : irrname(2)='Gam2'
420 0 : irrname(3)='Gam3'
421 0 : irrname(4)='Gam4'
422 0 : IF (soc) THEN
423 0 : nirr=8
424 0 : omega=EXP(CMPLX(0.0,-pi_const/4.0))
425 : char_table(5,1:4)=(/one, omega, -imi, &
426 0 : & -CONJG(omega)/)
427 : char_table(6,1:4)=(/one, CONJG(omega), imi, &
428 0 : & -omega/)
429 : char_table(7,1:4)=(/one, -omega, -imi, &
430 0 : & CONJG(omega)/)
431 : char_table(8,1:4)=(/one, -CONJG(omega), imi, &
432 0 : & omega/)
433 0 : irrname(5)='Gam5'
434 0 : irrname(6)='Gam6'
435 0 : irrname(7)='Gam7'
436 0 : irrname(8)='Gam8'
437 : ENDIF
438 0 : ELSE IF ((rot(2).EQ.1).AND.(rot(8).EQ.1)) THEN
439 0 : grpname='ZKUS'
440 : ENDIF
441 : CASE(5) ! D4, C4v, D2d, O and Td
442 0 : IF (rot(3).GT.0) THEN
443 0 : grpname='Td'
444 0 : char_table(1,1:5)=(/1.0, 1.0, 1.0, 1.0, 1.0/)
445 0 : char_table(2,1:5)=(/1.0, -1.0, 1.0, -1.0, 1.0/)
446 0 : char_table(3,1:5)=(/2.0, 0.0, 2.0, 0.0, -1.0/)
447 0 : char_table(4,1:5)=(/3.0, 1.0, -1.0, -1.0, 0.0/)
448 0 : char_table(5,1:5)=(/3.0, -1.0, -1.0, 1.0, 0.0/)
449 0 : irrname(1)='P1'
450 0 : irrname(2)='P2'
451 0 : irrname(3)='P3'
452 0 : irrname(4)='P5'
453 0 : irrname(5)='P4'
454 0 : ELSE IF (rot(10).GT.0) THEN
455 0 : grpname='D2d'
456 0 : char_table(1,1:5)=(/1.0, 1.0, 1.0, 1.0, 1.0/)
457 0 : char_table(2,1:5)=(/1.0, 1.0, -1.0, 1.0, -1.0/)
458 0 : char_table(3,1:5)=(/1.0, 1.0, 1.0, -1.0, -1.0/)
459 0 : char_table(4,1:5)=(/1.0, 1.0, -1.0, -1.0, 1.0/)
460 0 : char_table(5,1:5)=(/2.0, -2.0, 0.0, 0.0, 0.0/)
461 : ! standard order: E, 2, 2_x, -4, m
462 0 : irrname(1)='W1'
463 0 : irrname(2)='W2'
464 0 : irrname(3)='W1`'
465 0 : irrname(4)='W2`'
466 0 : irrname(5)='W3'
467 : ELSE
468 0 : grpname='C4v'
469 0 : char_table(1,1:5)=(/1.0, 1.0, 1.0, 1.0, 1.0/)
470 0 : char_table(2,1:5)=(/1.0, 1.0, 1.0, -1.0, -1.0/)
471 0 : char_table(3,1:5)=(/1.0, -1.0, 1.0, 1.0, -1.0/)
472 0 : char_table(4,1:5)=(/1.0, -1.0, 1.0, -1.0, 1.0/)
473 0 : char_table(5,1:5)=(/2.0, 0.0, -2.0, 0.0, 0.0/)
474 0 : irrname(1)='Del1'
475 0 : irrname(2)='Del1`'
476 0 : irrname(3)='Del2'
477 0 : irrname(4)='Del2`'
478 0 : irrname(5)='Del5'
479 : ENDIF
480 : CASE(6) ! C3i, D3d, C6, C3h, D6, C6v and D3h
481 0 : IF (rot(2).EQ.0) THEN
482 0 : IF (rot(7).EQ.0) THEN
483 0 : grpname='C3h'
484 0 : char_table(:,:) = 0.0
485 0 : WRITE(24,*) 'C3h: Character table missing in grp_k.F'
486 : ELSE
487 0 : grpname='C3i'
488 0 : omega=CMPLX(-0.5,SQRT(3./4.))
489 0 : char_table(1,1:6)=(/one, one, one, one, one, one/)
490 0 : char_table(2,1:6)=(/one, one, one,-one,-one,-one/)
491 : char_table(3,1:6)=(/one, omega, CONJG(omega),&
492 0 : & one, omega, CONJG(omega)/)
493 : char_table(4,1:6)=(/one, omega, CONJG(omega),&
494 0 : & -one,-omega,-CONJG(omega)/)
495 : char_table(5,1:6)=(/one, CONJG(omega), omega,&
496 0 : & one, CONJG(omega), omega/)
497 : char_table(6,1:6)=(/one, CONJG(omega), omega,&
498 0 : & -one,-CONJG(omega),-omega/)
499 0 : irrname(1)='Ag '
500 0 : irrname(2)='Au '
501 0 : irrname(3)='E1g '
502 0 : irrname(4)='E1u '
503 0 : irrname(5)='E2g '
504 0 : irrname(6)='E2u '
505 0 : IF (soc) THEN
506 0 : char_table(7,1:6)=(/one,-one, one, one,-one, one/)
507 0 : char_table(8,1:6)=(/one,-one, one,-one, one,-one/)
508 : char_table( 9,1:6)=(/one,-omega, CONJG(omega),&
509 0 : & one,-omega, CONJG(omega)/)
510 : char_table(10,1:6)=(/one,-omega, CONJG(omega),&
511 0 : & -one, omega,-CONJG(omega)/)
512 : char_table(11,1:6)=(/one,-CONJG(omega), omega,&
513 0 : & one,-CONJG(omega), omega/)
514 : char_table(12,1:6)=(/one,-CONJG(omega), omega,&
515 0 : & -one, CONJG(omega),-omega/)
516 0 : irrname(7)="Ag' "
517 0 : irrname(8)="Au' "
518 0 : irrname(9)="E1g'"
519 0 : irrname(10)="E1u'"
520 0 : irrname(11)="E2g'"
521 0 : irrname(12)="E2u'"
522 : ENDIF
523 : ENDIF
524 0 : ELSEIF (rot(2).GT.1) THEN
525 0 : grpname='D6'
526 0 : char_table(1,1:6)=(/1.0, 1.0, 1.0, 1.0, 1.0, 1.0/)
527 0 : char_table(2,1:6)=(/1.0, 1.0, 1.0, 1.0, -1.0, -1.0/)
528 0 : char_table(3,1:6)=(/1.0, -1.0, 1.0, -1.0, 1.0, -1.0/)
529 0 : char_table(4,1:6)=(/1.0, -1.0, 1.0, -1.0, -1.0, 1.0/)
530 0 : char_table(5,1:6)=(/2.0, 1.0, -1.0, -2.0, 0.0, 0.0/)
531 0 : char_table(6,1:6)=(/2.0, -1.0, -1.0, 2.0, 0.0, 0.0/)
532 0 : irrname(1)='Gam1'
533 0 : irrname(2)='Gam2'
534 0 : irrname(3)='Gam3'
535 0 : irrname(4)='Gam4'
536 0 : irrname(5)='Gam6'
537 0 : irrname(6)='Gam5'
538 : ELSE
539 0 : IF (rot(6).EQ.0) THEN
540 0 : grpname='D3d'
541 0 : char_table(1,1:6)=(/1.0, 1.0, 1.0, 1.0, 1.0, 1.0/)
542 0 : char_table(2,1:6)=(/1.0, 1.0, -1.0, 1.0, 1.0, -1.0/)
543 0 : char_table(3,1:6)=(/2.0, -1.0, 0.0, 2.0, -1.0, 0.0/)
544 0 : char_table(4,1:6)=(/1.0, 1.0, 1.0, -1.0, -1.0, -1.0/)
545 0 : char_table(5,1:6)=(/1.0, 1.0, -1.0, -1.0, -1.0, 1.0/)
546 0 : char_table(6,1:6)=(/2.0, -1.0, 0.0, -2.0, 1.0, 0.0/)
547 0 : irrname(1)='A1g'
548 0 : irrname(2)='A2g'
549 0 : irrname(3)='Eg'
550 0 : irrname(4)='A1u'
551 0 : irrname(5)='A2u'
552 0 : irrname(6)='Eu'
553 0 : ELSEIF (rot(6).GT.1) THEN
554 0 : grpname='C6'
555 0 : char_table(:,:) = 0.0
556 0 : WRITE(24,*) 'C6: Character table missing in grp_k.F'
557 : ELSE
558 0 : IF (members(3).GT.1) grpname='D3h' ! maybe this works
559 0 : IF (members(3).EQ.1) grpname='C6v'
560 0 : char_table(1,1:6)=(/1.0, 1.0, 1.0, 1.0, 1.0, 1.0/)
561 0 : char_table(2,1:6)=(/1.0, 1.0, 1.0, 1.0, -1.0, -1.0/)
562 0 : char_table(3,1:6)=(/1.0, -1.0, 1.0, -1.0, 1.0, -1.0/)
563 0 : char_table(4,1:6)=(/1.0, -1.0, 1.0, -1.0, -1.0, 1.0/)
564 0 : char_table(5,1:6)=(/2.0, 1.0, -1.0, -2.0, 0.0, 0.0/)
565 0 : char_table(6,1:6)=(/2.0, -1.0, -1.0, 2.0, 0.0, 0.0/)
566 0 : irrname(1)='Gam1'
567 0 : irrname(2)='Gam2'
568 0 : irrname(3)='Gam3'
569 0 : irrname(4)='Gam4'
570 0 : irrname(5)='Gam6'
571 0 : irrname(6)='Gam5'
572 : ENDIF
573 : ENDIF
574 : CASE(8) ! Th, C4h and D2h
575 0 : IF (rot(3).GT.0) THEN
576 0 : grpname='Th'
577 0 : char_table(:,:) = 0.0
578 0 : WRITE(24,*) 'Th: Character table missing in grp_k.F'
579 0 : ELSE IF (rot(2).EQ.3) THEN
580 0 : grpname='D2h'
581 0 : char_table(:,:) = 0.0
582 0 : WRITE(24,*) 'D2h: Character table missing in grp_k.F'
583 : ELSE
584 0 : grpname='C4h'
585 0 : char_table(:,:) = 0.0
586 0 : WRITE(24,*) 'C4h: Character table missing in grp_k.F'
587 : ENDIF
588 : CASE(10) ! D4h and Oh
589 0 : IF (rot(3).GT.0) THEN
590 0 : grpname='Oh'
591 : char_table(1,1:10)= (/1.0, 1.0, 1.0, 1.0, 1.0,&
592 0 : & 1.0, 1.0, 1.0, 1.0, 1.0/)
593 : char_table(2,1:10)= (/1.0, -1.0, 1.0, 1.0, -1.0,&
594 0 : & 1.0, -1.0, 1.0, 1.0, -1.0/)
595 : char_table(3,1:10)= (/2.0, 0.0, -1.0, 2.0, 0.0,&
596 0 : & 2.0, 0.0, -1.0, 2.0, 0.0/)
597 : char_table(4,1:10)= (/3.0, 1.0, 0.0, -1.0, -1.0,&
598 0 : & 3.0, 1.0, 0.0, -1.0, -1.0/)
599 : char_table(5,1:10)= (/3.0, -1.0, 0.0, -1.0, 1.0,&
600 0 : & 3.0, -1.0, 0.0, -1.0, 1.0/)
601 : char_table(6,1:10)= (/1.0, 1.0, 1.0, 1.0, 1.0,&
602 0 : & -1.0, -1.0, -1.0, -1.0, -1.0/)
603 : char_table(7,1:10)= (/1.0, -1.0, 1.0, 1.0, -1.0,&
604 0 : & -1.0, 1.0, -1.0, -1.0, 1.0/)
605 : char_table(8,1:10)= (/2.0, 0.0, -1.0, 2.0, 0.0,&
606 0 : & -2.0, 0.0, 1.0, -2.0, 0.0/)
607 : char_table(9,1:10)= (/3.0, 1.0, 0.0, -1.0, -1.0,&
608 0 : & -3.0, -1.0, 0.0, 1.0, 1.0/)
609 : char_table(10,1:10)=(/3.0, -1.0, 0.0, -1.0, 1.0,&
610 0 : & -3.0, 1.0, 0.0, 1.0, -1.0/)
611 :
612 0 : irrname(1)='Gam1+'
613 0 : irrname(2)='Gam2+'
614 0 : irrname(3)='Gam3+'
615 0 : irrname(4)='Gam4+'
616 0 : irrname(5)='Gam5+'
617 0 : irrname(6)='Gam1-'
618 0 : irrname(7)='Gam2-'
619 0 : irrname(8)='Gam3-'
620 0 : irrname(9)='Gam4-'
621 0 : irrname(10)='Gam5-'
622 : ELSE
623 0 : grpname='D4h'
624 : char_table(1,1:10)= (/1.0, 1.0, 1.0, 1.0, 1.0,&
625 0 : & 1.0, 1.0, 1.0, 1.0, 1.0/)
626 : char_table(2,1:10)= (/1.0, 1.0, 1.0, -1.0, -1.0,&
627 0 : & 1.0, 1.0, 1.0, -1.0, -1.0/)
628 : char_table(3,1:10)= (/1.0, -1.0, 1.0, 1.0, -1.0,&
629 0 : & 1.0, -1.0, 1.0, 1.0, -1.0/)
630 : char_table(4,1:10)= (/1.0, -1.0, 1.0, -1.0, 1.0,&
631 0 : & 1.0, -1.0, 1.0, -1.0, 1.0/)
632 : char_table(5,1:10)= (/2.0, 0.0, -2.0, 0.0, 0.0,&
633 0 : & 2.0, 0.0, -2.0, 0.0, 0.0/)
634 : char_table(6,1:10)= (/1.0, 1.0, 1.0, 1.0, 1.0,&
635 0 : & -1.0, -1.0, -1.0, -1.0, -1.0/)
636 : char_table(7,1:10)= (/1.0, 1.0, 1.0, -1.0, -1.0,&
637 0 : & -1.0, -1.0, -1.0, 1.0, 1.0/)
638 : char_table(8,1:10)= (/1.0, -1.0, 1.0, 1.0, -1.0,&
639 0 : & -1.0, 1.0, -1.0, -1.0, 1.0/)
640 : char_table(9,1:10)= (/1.0, -1.0, 1.0, -1.0, 1.0,&
641 0 : & -1.0, 1.0, -1.0, 1.0, -1.0/)
642 : char_table(10,1:10)=(/2.0, 0.0, -2.0, 0.0, 0.0,&
643 0 : & -2.0, 0.0, 2.0, 0.0, 0.0/)
644 :
645 : ! standard char-table assumes the order: E, 4, 2, 2_h, 2_h' where 4 and 2
646 : ! rotate around the same axis and 2_h is perpendicular to this axis. Check:
647 :
648 0 : CALL cross( rax(1,2),rax(1,3),kt )
649 0 : IF (kt(1)**2+kt(2)**2+kt(3)**2 < eps) THEN
650 : ! all ok
651 : ELSE
652 0 : CALL cross( rax(1,2),rax(1,4),kt )
653 0 : IF (kt(1)**2+kt(2)**2+kt(3)**2 < eps) THEN ! change 3 & 4
654 0 : CALL change_column(char_table,10,10,3,4)
655 : ! check, whether element 3 is perpendicular to element 2
656 0 : IF (DOT_PRODUCT(rax(1:3,2),rax(1:3,3)) < eps) THEN
657 : ! all ok
658 : ELSE ! change 4 & 5
659 0 : WRITE(*,*) DOT_PRODUCT(rax(1:3,2),rax(1:3,3))
660 0 : CALL change_column(char_table,10,10,4,5)
661 : ENDIF
662 : ELSE
663 0 : CALL cross( rax(1,2),rax(1,5),kt )
664 0 : IF (kt(1)**2+kt(2)**2+kt(3)**2 < eps) THEN ! change 3 & 5
665 0 : CALL change_column(char_table,10,10,3,5)
666 0 : IF (DOT_PRODUCT(rax(1:3,2),rax(1:3,3)) < eps) THEN ! change 4 & 5
667 0 : CALL change_column(char_table,10,10,4,5)
668 : ELSE
669 : ! all ok
670 : ENDIF
671 : ELSE
672 0 : CALL juDFT_error("D4h",calledby="grp_k")
673 : ENDIF
674 : ENDIF
675 : ENDIF
676 0 : irrname(1)='A1g'
677 0 : irrname(2)='A2g'
678 0 : irrname(3)='B1g'
679 0 : irrname(4)='B2g'
680 0 : irrname(5)='Eg'
681 0 : irrname(6)='A1u'
682 0 : irrname(7)='A2u'
683 0 : irrname(8)='B1u'
684 0 : irrname(9)='B2u'
685 0 : irrname(10)='Eu'
686 : ENDIF
687 : CASE(12) ! D6h and C6h
688 0 : IF (rot(2).EQ.3) THEN
689 0 : grpname='D6h'
690 0 : char_table(:,:) = 0.0
691 0 : WRITE(24,*) 'D6h: Character table missing in grp_k.F'
692 : ELSE
693 0 : grpname='C6h'
694 0 : char_table(:,:) = 0.0
695 0 : WRITE(24,*) 'C6h: Character table missing in grp_k.F'
696 : ENDIF
697 : CASE DEFAULT
698 0 : WRITE(24,*) 'Group of k not identified'
699 : END SELECT
700 :
701 : !>
702 :
703 0 : END SUBROUTINE grp_k
704 : !-------------------------------------------------------------------------------
705 :
706 : !<--SUBROUTINE rotaxis(alpha,beta,gamma,rax,theta)
707 0 : SUBROUTINE rotaxis(alpha,beta,gamma,rax,theta)
708 :
709 : ! Determines the rotation axis based on the Euler angles
710 : IMPLICIT NONE
711 : REAL, INTENT(IN) :: alpha,beta,gamma
712 : REAL, INTENT(OUT) :: rax(4),theta
713 :
714 : REAL :: costhe2,the2
715 : INTEGER :: i
716 :
717 0 : costhe2=COS(beta/2.0)*COS((alpha+gamma)/2.0)
718 0 : the2=ACOS(costhe2)
719 0 : rax=0.0
720 0 : IF (the2.GT.0.00001) THEN
721 0 : rax(1)=-SIN(beta/2.0)*SIN((alpha-gamma)/2.0)/SIN(the2)
722 0 : rax(2)= SIN(beta/2.0)*COS((alpha-gamma)/2.0)/SIN(the2)
723 0 : rax(3)= COS(beta/2.0)*SIN((alpha+gamma)/2.0)/SIN(the2)
724 : ENDIF
725 0 : DO i=1,3
726 0 : IF (ABS(rax(i)).GT.0.0001) THEN
727 0 : rax(i)=rax(i)/ABS(rax(i))
728 : ENDIF
729 : ENDDO
730 0 : rax(4)=rax(1)**2+rax(2)**2+rax(3)**2
731 0 : theta=the2*2.0
732 0 : END SUBROUTINE rotaxis
733 : !>
734 :
735 : !<--Subroutine euler
736 0 : SUBROUTINE euler(n,sym,cell,alpha,beta,gamma)
737 :
738 : ! determines the Euler angles corresponding the proper rotation part of mrot
739 :
740 : USE m_constants, ONLY : pi_const
741 : USE m_inv3
742 : USE m_types
743 : IMPLICIT NONE
744 : INTEGER,INTENT(IN) :: n
745 : TYPE(t_sym),INTENT(IN) :: sym
746 : TYPE(t_cell),INTENT(IN) :: cell
747 :
748 : REAL, INTENT(OUT) :: alpha,beta,gamma
749 :
750 : INTEGER :: det
751 : REAL :: mprop(3,3)
752 : REAL :: sina,sinb,sinc,cosa,cosb,cosc
753 : REAL :: amatinv(3,3),detr
754 :
755 0 : CALL inv3(cell%amat,amatinv,detr)
756 : det= sym%mrot(1,1,n)*sym%mrot(2,2,n)*sym%mrot(3,3,n) +&
757 : & sym%mrot(1,2,n)*sym%mrot(2,3,n)*sym%mrot(3,1,n) +&
758 : & sym%mrot(2,1,n)*sym%mrot(3,2,n)*sym%mrot(1,3,n) -&
759 : & sym%mrot(1,3,n)*sym%mrot(2,2,n)*sym%mrot(3,1,n) -&
760 : & sym%mrot(2,3,n)*sym%mrot(3,2,n)*sym%mrot(1,1,n) -&
761 0 : & sym%mrot(2,1,n)*sym%mrot(1,2,n)*sym%mrot(3,3,n)
762 :
763 : ! Take the proper rotation
764 0 : mprop=REAL(det)*MATMUL(cell%amat,MATMUL(REAL(sym%mrot(:,:,n)),amatinv))
765 : ! Euler angles
766 0 : cosb = mprop(3,3)
767 0 : sinb = 1.00 - cosb*cosb
768 0 : sinb = MAX(sinb,0.00)
769 0 : sinb = SQRT(sinb)
770 : !
771 : ! if beta = 0 or pi , only alpha+gamma or -gamma have a meaning:
772 : !
773 0 : IF ( ABS(sinb).LT.1.0e-5 ) THEN
774 0 : beta = 0.0
775 0 : IF ( cosb.LT.0.0 ) beta = pi_const
776 0 : gamma = 0.0
777 0 : cosa = mprop(1,1)/cosb
778 0 : sina = mprop(1,2)/cosb
779 0 : IF ( ABS(sina).LT.1.0e-5 ) THEN
780 0 : alpha=0.0
781 0 : IF ( cosa.LT.0.0 ) alpha=alpha+pi_const
782 : ELSE
783 0 : alpha = 0.5*pi_const - ATAN(cosa/sina)
784 0 : IF ( sina.LT.0.0 ) alpha=alpha+pi_const
785 : ENDIF
786 : ELSE
787 0 : beta = 0.5*pi_const - ATAN(cosb/sinb)
788 : !
789 : ! determine alpha and gamma from d13 d23 d32 d31
790 : !
791 0 : cosa = mprop(3,1)/sinb
792 0 : sina = mprop(3,2)/sinb
793 0 : cosc =-mprop(1,3)/sinb
794 0 : sinc = mprop(2,3)/sinb
795 0 : IF ( ABS(sina).LT.1.0e-5 ) THEN
796 0 : alpha=0.0
797 0 : IF ( cosa.LT.0.0 ) alpha=alpha+pi_const
798 : ELSE
799 0 : alpha = 0.5*pi_const - ATAN(cosa/sina)
800 0 : IF ( sina.LT.0.0 ) alpha=alpha+pi_const
801 : ENDIF
802 0 : IF ( ABS(sinc).LT.1.0e-5 ) THEN
803 0 : gamma = 0.0
804 0 : IF ( cosc.LT.0.0 ) gamma=gamma+pi_const
805 : ELSE
806 0 : gamma = 0.5*pi_const - ATAN(cosc/sinc)
807 0 : IF ( sinc.LT.0.0 ) gamma=gamma+pi_const
808 : ENDIF
809 :
810 : ENDIF
811 :
812 0 : END SUBROUTINE euler
813 : !>
814 : !-----------------------------------------------------------
815 0 : SUBROUTINE change_column( char_table,nx,ny,r1,r2 )
816 :
817 : INTEGER, INTENT (IN) :: nx,ny,r1,r2
818 : COMPLEX, INTENT(INOUT) :: char_table(:,:)
819 :
820 0 : COMPLEX tmpc(nx)
821 :
822 0 : tmpc(1:nx) = char_table(1:nx,r1)
823 0 : char_table(1:nx,r1) = char_table(1:nx,r2)
824 0 : char_table(1:nx,r2) = tmpc(1:nx)
825 0 : tmpc(1:nx) = char_table(1:nx,r1+nx/2)
826 0 : char_table(1:nx,r1+nx/2) = char_table(1:nx,r2+nx/2)
827 0 : char_table(1:nx,r2+nx/2) = tmpc(1:nx)
828 :
829 0 : END SUBROUTINE change_column
830 0 : END MODULE m_grp_k
|
