Line data Source code
1 : !!****m* ABINIT/m_ptgroups
2 : !! NAME
3 : !! m_ptgroups
4 : !!
5 : !! FUNCTION
6 : !! This module contains the irreducible representations and the
7 : !! character tables of the 32 point groups.
8 : !!
9 : !! COPYRIGHT
10 : !! Copyright (C) 2010-2026 ABINIT group (MG)
11 : !! This file is distributed under the terms of the
12 : !! GNU General Public License, see ~abinit/COPYING
13 : !! or http://www.gnu.org/copyleft/gpl.txt .
14 : !!
15 :
16 : #if defined HAVE_CONFIG_H
17 : #include "config.h"
18 : #endif
19 :
20 : #include "abi_common.h"
21 :
22 : module m_ptgroups
23 :
24 : use defs_basis
25 : use m_errors
26 : use m_abicore
27 : use m_defs_ptgroups
28 :
29 : use m_io_tools, only : open_file
30 : use m_fstrings, only : sjoin
31 : use m_numeric_tools, only : get_trace, cmplx_sphcart
32 : use m_matrix, only : mati3inv
33 : use m_symtk, only : sg_multable
34 :
35 : ! Import group tables
36 : use m_ptg_C1
37 : use m_ptg_Ci
38 : use m_ptg_C2
39 : use m_ptg_Cs
40 : use m_ptg_C2h
41 : use m_ptg_D2
42 : use m_ptg_C2v
43 : use m_ptg_D2h
44 : use m_ptg_C4
45 : use m_ptg_S4
46 : use m_ptg_C4h
47 : use m_ptg_D4
48 : use m_ptg_C4v
49 : use m_ptg_D2d
50 : use m_ptg_D4h
51 : use m_ptg_C3
52 : use m_ptg_C3i
53 : use m_ptg_D3
54 : use m_ptg_C3v
55 : use m_ptg_D3d
56 : use m_ptg_C6
57 : use m_ptg_C3h
58 : use m_ptg_C6h
59 : use m_ptg_D6
60 : use m_ptg_C6v
61 : use m_ptg_D3h
62 : use m_ptg_D6h
63 : use m_ptg_T
64 : use m_ptg_Th
65 : use m_ptg_O
66 : use m_ptg_Td
67 : use m_ptg_Oh
68 :
69 : implicit none
70 :
71 : private
72 :
73 : public :: get_point_group
74 : public :: get_classes
75 : public :: show_character_tables
76 : public :: point_group_init
77 : public :: sum_irreps
78 : public :: mult_table
79 : !public :: polish_irreps
80 :
81 :
82 : contains
83 : !!***
84 :
85 : !----------------------------------------------------------------------
86 :
87 : !!****f* m_ptgroups/get_point_group
88 : !! NAME
89 : !! get_point_group
90 : !!
91 : !! FUNCTION
92 : !!
93 : !! INPUTS
94 : !! ptg_name=point group name as returned by symptgroup.
95 : !!
96 : !! OUTPUT
97 : !! nsym=Number of symmetries in the point group.
98 : !! nclass=Number of classes.
99 : !! sym(3,3,nsym)=Elements of the point group ordered by classe.
100 : !! class_ids(2,nclass)=Initial and final index in sym, for each
101 : !! Irreps(nclass)=Datatype gathering data on the different irreducible representations.
102 : !!
103 : !! SOURCE
104 :
105 0 : subroutine get_point_group(ptg_name, nsym, nclass, sym, class_ids, class_names, Irreps)
106 :
107 : !Arguments ------------------------------------
108 : !scalars
109 : character(len=*),intent(in) :: ptg_name
110 : integer,intent(out) :: nsym, nclass
111 : !arrays
112 : integer,allocatable,intent(out) :: sym(:,:,:),class_ids(:,:)
113 : character(len=5),allocatable,intent(out) :: class_names(:)
114 : type(irrep_t),allocatable,intent(out) :: Irreps(:)
115 :
116 : !Local variables-------------------------------
117 : integer :: irp,isym
118 : ! *************************************************************************
119 :
120 0 : SELECT CASE (TRIM(ADJUSTL(ptg_name)))
121 : CASE ('1')
122 0 : call ptg_C1 (nsym,nclass,sym,class_ids,class_names,Irreps)
123 : CASE ('-1')
124 0 : call ptg_Ci (nsym,nclass,sym,class_ids,class_names,Irreps)
125 : CASE ('2')
126 0 : call ptg_C2 (nsym,nclass,sym,class_ids,class_names,Irreps)
127 : CASE ('m',"-2") ! Abinit uses "-2"
128 0 : call ptg_Cs (nsym,nclass,sym,class_ids,class_names,Irreps)
129 : CASE ('2/m')
130 0 : call ptg_C2h (nsym,nclass,sym,class_ids,class_names,Irreps)
131 : CASE ('222')
132 0 : call ptg_D2 (nsym,nclass,sym,class_ids,class_names,Irreps)
133 : CASE ('mm2')
134 0 : call ptg_C2v (nsym,nclass,sym,class_ids,class_names,Irreps)
135 : CASE ('mmm')
136 0 : call ptg_D2h (nsym,nclass,sym,class_ids,class_names,Irreps)
137 : CASE ('4')
138 0 : call ptg_C4 (nsym,nclass,sym,class_ids,class_names,Irreps)
139 : CASE ('-4')
140 0 : call ptg_S4 (nsym,nclass,sym,class_ids,class_names,Irreps)
141 : CASE ('4/m')
142 0 : call ptg_C4h (nsym,nclass,sym,class_ids,class_names,Irreps)
143 : CASE ('422')
144 0 : call ptg_D4 (nsym,nclass,sym,class_ids,class_names,Irreps)
145 : CASE ('4mm')
146 0 : call ptg_C4v (nsym,nclass,sym,class_ids,class_names,Irreps)
147 : CASE ('-42m')
148 0 : call ptg_D2d (nsym,nclass,sym,class_ids,class_names,Irreps)
149 : CASE ('4/mmm')
150 0 : call ptg_D4h (nsym,nclass,sym,class_ids,class_names,Irreps)
151 : CASE ('3')
152 0 : call ptg_C3 (nsym,nclass,sym,class_ids,class_names,Irreps)
153 : CASE ('-3')
154 0 : call ptg_C3i (nsym,nclass,sym,class_ids,class_names,Irreps)
155 : CASE ('32')
156 0 : call ptg_D3 (nsym,nclass,sym,class_ids,class_names,Irreps)
157 : CASE ('3m')
158 0 : call ptg_C3v (nsym,nclass,sym,class_ids,class_names,Irreps)
159 : CASE ('-3m')
160 0 : call ptg_D3d (nsym,nclass,sym,class_ids,class_names,Irreps)
161 : CASE ('6')
162 0 : call ptg_C6 (nsym,nclass,sym,class_ids,class_names,Irreps)
163 : CASE ('-6')
164 0 : call ptg_C3h (nsym,nclass,sym,class_ids,class_names,Irreps)
165 : CASE ('6/m')
166 0 : call ptg_C6h (nsym,nclass,sym,class_ids,class_names,Irreps)
167 : CASE ('622')
168 0 : call ptg_D6 (nsym,nclass,sym,class_ids,class_names,Irreps)
169 : CASE ('6mm')
170 0 : call ptg_C6v (nsym,nclass,sym,class_ids,class_names,Irreps)
171 : CASE ('-62m')
172 0 : call ptg_D3h (nsym,nclass,sym,class_ids,class_names,Irreps)
173 : CASE ('6/mmm')
174 0 : call ptg_D6h (nsym,nclass,sym,class_ids,class_names,Irreps)
175 : CASE ('23')
176 0 : call ptg_T (nsym,nclass,sym,class_ids,class_names,Irreps)
177 : CASE ('m-3')
178 0 : call ptg_Th (nsym,nclass,sym,class_ids,class_names,Irreps)
179 : CASE ('432')
180 0 : call ptg_O (nsym,nclass,sym,class_ids,class_names,Irreps)
181 : CASE ('-43m')
182 0 : call ptg_Td (nsym,nclass,sym,class_ids,class_names,Irreps)
183 : CASE ('m-3m')
184 0 : call ptg_Oh (nsym,nclass,sym,class_ids,class_names,Irreps)
185 : CASE DEFAULT
186 0 : ABI_BUG(sjoin("Unknown value for ptg_name:", ptg_name))
187 : END SELECT
188 :
189 : ! Calculate the trace of each irreducible representation in order to have the character at hand.
190 0 : do irp=1,SIZE(Irreps)
191 0 : ABI_MALLOC(Irreps(irp)%trace, (nsym))
192 0 : do isym=1,nsym
193 0 : Irreps(irp)%trace(isym) = get_trace(Irreps(irp)%mat(:,:,isym))
194 : end do
195 : end do
196 :
197 0 : end subroutine get_point_group
198 : !!***
199 :
200 : !----------------------------------------------------------------------
201 :
202 : !!****f* m_ptgroups/get_classes
203 : !! NAME
204 : !! get_classes
205 : !!
206 : !! FUNCTION
207 : !! Given a set of nsym 3x3 operations in reciprocal or real space,
208 : !! which are supposed to form a group, this routine divides the group into classes.
209 : !!
210 : !! INPUTS
211 : !! nsym=number of symmetry operation.
212 : !! sym(3,3,nsym)=The symmetry operations.
213 : !!
214 : !! OUTPUT
215 : !! nclass=The number of classes
216 : !! nelements(1:nclass)=For each class, the number of elements
217 : !! elements_idx(ii,1:nclass)=For each class, this table gives the index of its elements (ii=1,..,nelements(iclass))
218 : !!
219 : !! NOTES
220 : !! * A class is defined as the set of distinct elements obtained by
221 : !! considering for each element, S, of the group all its conjugate
222 : !! elements X^-1 S X where X ranges over all the elements of the group.
223 : !!
224 : !! * It does not work in case of non-collinear magnetism.
225 : !!
226 : !! * The routine assumes that anti-ferromagnetic symmetries (if any) have been removed by the caller.
227 : !!
228 : !! SOURCE
229 :
230 0 : subroutine get_classes(nsym, sym, nclass, nelements, elements_idx)
231 :
232 : !Arguments ------------------------------------
233 : !scalars
234 : integer,intent(in) :: nsym
235 : integer,intent(out) :: nclass
236 : !arrays
237 : integer,intent(in) :: sym(3,3,nsym)
238 : integer,intent(out) :: nelements(nsym),elements_idx(nsym,nsym)
239 :
240 : !Local variables-------------------------------
241 : !scalars
242 : integer :: isym,jsym,ksym,identity_idx,ierr
243 : character(len=500) :: msg
244 : !arrays
245 0 : integer :: cjg(3,3),ss(3,3),xx(3,3),xxm1(3,3),test(3,3), identity(3,3), dummy_symafm(nsym)
246 0 : logical :: found(nsym),found_identity
247 : !************************************************************************
248 :
249 : ! Check if identity is present in the first position
250 0 : identity=RESHAPE((/1,0,0,0,1,0,0,0,1/),(/3,3/)); found_identity=.FALSE.
251 :
252 0 : do isym=1,nsym
253 0 : if (all(sym(:,:,isym) == identity)) then
254 0 : found_identity=.TRUE.; identity_idx=isym; exit
255 : end if
256 : end do
257 :
258 0 : if (.not. found_identity .or. identity_idx /= 1) then
259 : write(msg, '(3a)')&
260 0 : 'Either identity is not present or it is not the first operation ',ch10,'Check set of symmetry operations '
261 0 : ABI_ERROR(msg)
262 : end if
263 :
264 0 : dummy_symafm = 1
265 0 : call sg_multable(nsym, dummy_symafm, sym, ierr)
266 0 : ABI_CHECK(ierr == 0, "Error in group closure")
267 :
268 0 : nclass=0; nelements(:)=0; elements_idx(:,:)=0; found(:)=.FALSE.
269 0 : do isym=1,nsym
270 0 : if (.not.found(isym)) then
271 0 : nclass=nclass+1
272 0 : ss(:,:)=sym(:,:,isym)
273 :
274 0 : do jsym=1,nsym ! Form conjugate.
275 0 : xx(:,:)=sym(:,:,jsym)
276 0 : call mati3inv(xx,xxm1) ; xxm1=TRANSPOSE(xxm1)
277 0 : cjg(:,:)=MATMUL(xxm1,MATMUL(ss,xx))
278 0 : do ksym=1,nsym ! Is it already found?
279 0 : test(:,:)=sym(:,:,ksym)
280 0 : if (.not.found(ksym).and.(ALL((test-cjg)==0))) then
281 0 : found(ksym)=.TRUE.
282 0 : nelements(nclass)=nelements(nclass)+1
283 0 : elements_idx(nelements(nclass),nclass)=ksym
284 : end if
285 : end do
286 : end do
287 :
288 : end if
289 : end do
290 :
291 0 : end subroutine get_classes
292 : !!***
293 :
294 : !----------------------------------------------------------------------
295 :
296 : !!****f* m_ptgroups/show_character_tables
297 : !! NAME
298 : !! show_character_tables
299 : !!
300 : !! FUNCTION
301 : !! Printout of the caracter tables of the 32 point groups.
302 : !!
303 : !! INPUTS
304 : !! [unit]=Unit number of output file. Defaults to std_out
305 : !!
306 : !! OUTPUT
307 : !! Only writing.
308 : !!
309 : !! SOURCE
310 :
311 0 : subroutine show_character_tables(unit)
312 :
313 : !Arguments ------------------------------------
314 : !scalars
315 : integer,optional,intent(in) :: unit
316 :
317 : !Local variables-------------------------------
318 : integer :: igrp,my_unt
319 : character(len=5) :: ptg_name
320 0 : type(point_group_t) :: Ptg
321 : !arrays
322 : !integer,allocatable :: elements_idx(:,:),nelements(:)
323 : ! *********************************************************************
324 :
325 0 : my_unt = std_out; if (PRESENT(unit)) my_unt=unit
326 :
327 0 : do igrp=1,SIZE(ptgroup_names)
328 0 : ptg_name = ptgroup_names(igrp)
329 0 : call point_group_init(Ptg, ptg_name)
330 0 : call Ptg%print(unit=my_unt)
331 : !allocate(nelements(Ptg%nsym),elements_idx(Ptg%nsym,Ptg%nsym))
332 : !call get_classes(Ptg%nsym,Ptg%sym,nclass,nelements,elements_idx)
333 : !deallocate(nelements,elements_idx)
334 0 : call Ptg%free()
335 : end do
336 :
337 0 : end subroutine show_character_tables
338 : !!***
339 :
340 : !----------------------------------------------------------------------
341 :
342 : !!****f* m_ptgroups/point_group_init
343 : !! NAME
344 : !! point_group_init
345 : !!
346 : !! FUNCTION
347 : !! Creation method for the point_group_t datatype.
348 : !!
349 : !! INPUTS
350 : !! ptg_name=The name of the point group (International conventions).
351 : !!
352 : !! OUTPUT
353 : !! The datatype completely initialized.
354 : !!
355 : !! SOURCE
356 :
357 0 : subroutine point_group_init(Ptg, ptg_name)
358 :
359 : !Arguments ------------------------------------
360 : class(point_group_t),intent(inout) :: Ptg
361 : character(len=5),intent(in) :: ptg_name
362 : ! *********************************************************************
363 :
364 0 : Ptg%gname = ptg_name
365 0 : call get_point_group(Ptg%gname,Ptg%nsym,Ptg%nclass,Ptg%sym,Ptg%class_ids,Ptg%class_names,Ptg%Irreps)
366 :
367 0 : end subroutine point_group_init
368 : !!***
369 :
370 : !----------------------------------------------------------------------
371 :
372 : !!****f* m_ptgroups/mult_table
373 : !! NAME
374 : !! mult_table
375 : !!
376 : !! FUNCTION
377 : !! Given a set of nsym 3x3 operations which are supposed to form a group,
378 : !! this routine constructs the multiplication table of the group.
379 : !!
380 : !! INPUTS
381 : !! nsym=number of symmetry operation
382 : !! sym(3,3,nsym)=the operations
383 : !!
384 : !! OUTPUT
385 : !! mtab(nsym,nsym)=The index of the product S_i * S_j in the input set sym.
386 : !!
387 : !! SOURCE
388 :
389 0 : subroutine mult_table(nsym, sym, mtab)
390 :
391 : !Arguments ------------------------------------
392 : !scalars
393 : integer,intent(in) :: nsym
394 : !arrays
395 : integer,intent(in) :: sym(3,3,nsym)
396 : integer,intent(out) :: mtab(nsym,nsym)
397 :
398 : !Local variables-------------------------------
399 : !scalars
400 : integer :: isym,jsym,ksym
401 : !character(len=500) :: msg
402 : !arrays
403 0 : integer :: prod_ij(3,3),found(nsym)
404 : !************************************************************************
405 :
406 0 : do jsym=1,nsym
407 0 : found(:)=0 ! Each symmetry should compare only once in a given (row|col).
408 :
409 0 : do isym=1,nsym
410 0 : prod_ij = MATMUL(sym(:,:,isym),sym(:,:,jsym))
411 0 : do ksym=1,nsym
412 0 : if ( ALL(prod_ij == sym(:,:,ksym)) ) then
413 0 : found(ksym)=found(ksym)+1
414 0 : mtab(isym,jsym) = ksym
415 : end if
416 : end do
417 : end do ! jsym
418 :
419 0 : if (any(found /= 1)) then
420 0 : write(std_out,*)"found = ",found
421 0 : ABI_ERROR("Input elements do not form a group")
422 : end if
423 : end do ! isym
424 :
425 0 : end subroutine mult_table
426 : !!***
427 :
428 : !----------------------------------------------------------------------
429 :
430 : !!****f* m_ptgroups/groupk_from_file
431 : !! NAME
432 : !! groupk_from_file
433 : !!
434 : !! FUNCTION
435 : !! Initialize the group_k_t datatype from an external database retrieved from
436 : !! the Bilbao server via the ptg.py script.
437 : !!
438 : !! INPUTS
439 : !! fname(len=*)=file name
440 : !!
441 : !! OUTPUT
442 : !! ierr=Status error
443 : !! Lgrps<group_k_t>=The structure completely initialized.
444 : !!
445 : !! TODO
446 : !! This is a stub. I still have to complete the fileformat for the Bilbao database.
447 : !!
448 : !! SOURCE
449 :
450 0 : subroutine groupk_from_file(Lgrps, spgroup, fname, nkpt, klist, ierr)
451 :
452 : !Arguments ------------------------------------
453 : !scalars
454 : integer,intent(in) :: spgroup
455 : integer,intent(out) :: ierr,nkpt
456 : character(len=fnlen),intent(in) :: fname
457 : !arrays
458 : type(group_k_t),target,allocatable :: Lgrps(:)
459 : real(dp),pointer :: klist(:,:)
460 :
461 : !Local variables-------------------------------
462 : !scalars
463 : integer,parameter :: last_file_version=1
464 : integer :: unt,fvers,ik,nsym_ltgk,isym,nirreps_k,irp,icls
465 : integer :: irrep_idx,irrep_dim,sym_idx,ita_spgroup,nels,prev,now
466 : character(len=IRREPNAME_LEN) :: irrep_name
467 : character(len=1) :: basis
468 : character(len=500) :: msg
469 : type(group_k_t),pointer :: Gk
470 : !arrays
471 0 : integer,allocatable :: nelements(:),elements_idx(:,:)
472 : real(dp) :: kpt(3)
473 0 : character(len=10),allocatable :: kname(:)
474 : type(irrep_t),pointer :: OneIrr
475 : ! *************************************************************************
476 :
477 0 : ierr=0
478 0 : if (open_file(fname,msg,newunit=unt,form="formatted") /=0) then
479 0 : ABI_ERROR(msg)
480 : end if
481 :
482 0 : read(unt,*,ERR=10) ! Skip the header.
483 0 : read(unt,*,ERR=10) fvers ! File version.
484 0 : if (fvers > last_file_version) then
485 0 : write(msg,"(2(a,i0))")" Found file format= ",fvers," but the latest supported version is: ",last_file_version
486 0 : ABI_ERROR(msg)
487 : end if
488 0 : read(unt,*,ERR=10) ita_spgroup
489 0 : read(unt,*,ERR=10) basis
490 :
491 0 : if (spgroup/=ita_spgroup) then
492 0 : write(msg,'(a,2i0)')" Input space group does not match with the value reported on file: ",spgroup,ita_spgroup
493 0 : ABI_ERROR(msg)
494 : end if
495 :
496 0 : if (basis /= "b") then
497 0 : ABI_ERROR(" Wrong value for basis: "//TRIM(basis))
498 : end if
499 :
500 : ! Read the list of the k-points.
501 0 : read(unt,*,ERR=10) nkpt
502 :
503 0 : ABI_MALLOC(Lgrps,(nkpt))
504 0 : ABI_MALLOC(klist,(3,nkpt))
505 0 : ABI_MALLOC(kname,(nkpt))
506 0 : do ik=1,nkpt
507 0 : read(unt,*,ERR=10) klist(:,ik), kname(ik)
508 : end do
509 :
510 : ! Read tables for each k-point
511 0 : do ik=1,nkpt
512 :
513 0 : read(unt,*,ERR=10) kpt
514 0 : read(unt,*,ERR=10) nsym_ltgk
515 0 : Gk => Lgrps(ik)
516 :
517 0 : Gk%spgroup = ita_spgroup
518 0 : Gk%nsym = nsym_ltgk
519 0 : Gk%point = kpt
520 0 : ABI_MALLOC(Gk%sym,(3,3,nsym_ltgk))
521 0 : ABI_MALLOC(Gk%tnons,(3,nsym_ltgk))
522 :
523 0 : do isym=1,nsym_ltgk ! Read symmetries of the little group.
524 0 : read(unt,*,ERR=10) Gk%sym(:,:,isym)
525 0 : read(unt,*,ERR=10) Gk%tnons(:,isym)
526 : end do
527 :
528 0 : ABI_MALLOC(nelements,(nsym_ltgk))
529 0 : ABI_MALLOC(elements_idx,(nsym_ltgk,nsym_ltgk))
530 :
531 0 : call get_classes(nsym_ltgk,Gk%sym,Gk%nclass,nelements,elements_idx)
532 :
533 : ! The operations reported on the file are supposed to be packed in classes
534 : ! otherwise one should perform a rearrangement of the indices.
535 0 : prev = 0
536 0 : do icls=1,Gk%nclass
537 0 : do isym=1,nelements(icls)
538 0 : now = elements_idx(isym,icls)
539 0 : if ( (now-prev) /= 1 ) then
540 0 : write(msg,"(2(a,i0))")" Symmetries on file are not ordered in classes. icls= ",icls,", isym= ",isym
541 0 : ABI_ERROR(msg)
542 : else
543 0 : prev = now
544 : end if
545 : end do
546 : end do
547 :
548 0 : ABI_MALLOC(Gk%class_ids,(2,Gk%nclass))
549 0 : do icls=1,Gk%nclass
550 0 : nels = nelements(icls)
551 0 : Gk%class_ids(1,icls) = elements_idx(1, icls)
552 0 : Gk%class_ids(2,icls) = elements_idx(nels,icls)
553 : end do
554 :
555 0 : ABI_FREE(nelements)
556 0 : ABI_FREE(elements_idx)
557 :
558 : ! Read the irreducible representations.
559 0 : read(unt,*,ERR=10) nirreps_k
560 0 : ABI_CHECK(Gk%nclass == nirreps_k,"Gk%nclass /= nirreps_k")
561 :
562 : !$$ allocate(Gk%class_names(Gk%nclass))
563 0 : ABI_MALLOC(Gk%Irreps, (nirreps_k))
564 :
565 0 : do irp=1,nirreps_k
566 0 : OneIrr => Gk%Irreps(irp)
567 0 : read(unt,*,ERR=10) irrep_idx, irrep_dim, irrep_name
568 0 : call OneIrr%init(nsym_ltgk, irrep_dim, irrep_name)
569 0 : do isym=1,nsym_ltgk
570 0 : read(unt,*,ERR=10) sym_idx, OneIrr%mat(:,:,isym)
571 0 : ABI_CHECK(sym_idx==irp,"sym_idx/=irp!")
572 : ! Matrix elements on file are in the form (rho, theta) with theta given in degrees.
573 0 : call cmplx_sphcart(OneIrr%mat(:,:,isym),from="Sphere",units="Degrees")
574 0 : OneIrr%trace(isym) = get_trace(OneIrr%mat(:,:,isym))
575 : end do
576 : end do
577 :
578 : end do
579 :
580 0 : close(unt)
581 0 : RETURN
582 :
583 : ! Handle IO-error.
584 0 : 10 ierr=1
585 0 : close(unt)
586 0 : RETURN
587 :
588 0 : end subroutine groupk_from_file
589 : !!***
590 :
591 : !----------------------------------------------------------------------
592 :
593 : !!****f* m_ptgroups/sum_irreps
594 : !! NAME
595 : !! sum_irreps
596 : !!
597 : !! FUNCTION
598 : !!
599 : !! INPUTS
600 : !!
601 : !! OUTPUT
602 : !!
603 : !! SOURCE
604 :
605 0 : function sum_irreps(Irrep1,Irrep2,ii,jj,kk,ll) result(res)
606 :
607 : !Arguments ------------------------------------
608 : !scalars
609 : class(irrep_t),intent(in) :: Irrep1,Irrep2
610 : integer,intent(in) :: ii,jj,kk,ll
611 : !arrays
612 : complex(dp) :: res
613 :
614 : !Local variables-------------------------------
615 : integer :: isym,nsym,ierr
616 : !character(len=500) :: msg
617 : ! *********************************************************************
618 :
619 0 : ierr=0; res = czero
620 :
621 0 : nsym = Irrep1%nsym
622 0 : if (nsym /= Irrep2%nsym) then
623 0 : ABI_WARNING("Irreps have different nsym")
624 0 : ierr=ierr+1
625 : end if
626 :
627 0 : if (Irrep1%dim /= Irrep2%dim) then
628 0 : ABI_WARNING("Irreps have different dimensions")
629 0 : write(std_out,*)Irrep1%dim,Irrep2%dim
630 0 : ierr=ierr+1
631 : end if
632 :
633 0 : if (ii > Irrep2%dim .or. jj > Irrep2%dim .or. kk > Irrep1%dim .or. ll > Irrep1%dim) then
634 0 : ABI_WARNING("Wrong indices")
635 0 : write(std_out,*)ii,Irrep2%dim,jj,Irrep2%dim,kk>Irrep1%dim,ll,Irrep1%dim
636 0 : ierr=ierr+1
637 : end if
638 :
639 0 : if (ierr/=0) RETURN
640 :
641 0 : do isym=1,nsym
642 0 : res = res + DCONJG(Irrep1%mat(ii,jj,isym)) * Irrep2%mat(kk,ll,isym)
643 : end do
644 :
645 0 : end function sum_irreps
646 : !!***
647 :
648 : end module m_ptgroups
649 : !!***
|