TABLE OF CONTENTS


ABINIT/m_ptgroups [ Modules ]

[ Top ] [ Modules ]

NAME

 m_ptgroups

FUNCTION

  This module contains the irreducible representations and the
  character tables of the 32 point groups.

COPYRIGHT

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


m_ptgroups/copy_irrep [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

  copy_irrep

FUNCTION

  Perform a copy of a set of irrep_t datatypes. Optionally one can multiply
  by a phase factor.

PARENTS

      m_esymm

CHILDREN

      irrep_free

SOURCE

1038 subroutine copy_irrep(In_irreps,Out_irreps,phase_fact)
1039 
1040 
1041 !This section has been created automatically by the script Abilint (TD).
1042 !Do not modify the following lines by hand.
1043 #undef ABI_FUNC
1044 #define ABI_FUNC 'copy_irrep'
1045 !End of the abilint section
1046 
1047  implicit none
1048 
1049 !Arguments ------------------------------------
1050  type(irrep_t),intent(in) :: In_irreps(:)
1051  type(irrep_t),intent(inout) :: Out_irreps(:)
1052  complex(dpc),optional,intent(in) :: phase_fact(:)
1053 
1054 !Local variables-------------------------------
1055 !scalars
1056  integer :: irp,dim1,dim2,in_nsym,in_dim,isym
1057  character(len=500) :: msg
1058 !arrays
1059  complex(dpc) :: my_phase_fact(In_irreps(1)%nsym)
1060 ! *********************************************************************
1061 
1062  !@irrep_t
1063  dim1 = SIZE( In_irreps)
1064  dim2 = SIZE(Out_irreps)
1065  if (dim1/=dim2) then
1066    msg = " irreps to be copied have different dimension"
1067    MSG_ERROR(msg)
1068  end if
1069 
1070  my_phase_fact=cone
1071  if (PRESENT(phase_fact)) then
1072    my_phase_fact=phase_fact
1073    if (SIZE(phase_fact) /= In_irreps(1)%nsym) then
1074      msg = " irreps to be copied have different dimension"
1075      MSG_ERROR(msg)
1076    end if
1077  end if
1078 
1079  do irp=1,dim1
1080    in_dim  = In_irreps(irp)%dim
1081    in_nsym = In_irreps(irp)%nsym
1082    call init_irrep(Out_irreps(irp),in_nsym,in_dim)
1083    Out_irreps(irp)%name = In_irreps(irp)%name
1084    do isym=1,in_nsym
1085      Out_irreps(irp)%mat(:,:,isym) = In_irreps(irp)%mat(:,:,isym) * my_phase_fact(isym)
1086      Out_irreps(irp)%trace(isym) = get_trace(Out_irreps(irp)%mat(:,:,isym))
1087    end do
1088  end do
1089 
1090 end subroutine copy_irrep

m_ptgroups/get_classes [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

 get_classes

FUNCTION

  Given a set of nsym 3x3 operations in reciprocal or real space,
  which are supposed to form a group, this routine divides the group into classes.

INPUTS

 nsym=number of symmetry operation.
 sym(3,3,nsym)=The symmetry operations.

OUTPUT

 nclass=The number of classes
 nelements(1:nclass)=For each class, the number of elements
 elements_idx(ii,1:nclass)=For each class, this table gives the index
   of its elements (ii=1,..,nelements(iclass))

NOTES

  * A class is defined as the set of distinct elements obtained by
    considering for each element, S, of the group all its conjugate
    elements X^-1 S X where X range over all the elements of the group.

  * It does not work in case of non-collinear magnetism.

  * The routine assumes that anti-ferromagnetic symmetries (if any) have been removed by the caller.

PARENTS

      m_esymm,m_ptgroups

CHILDREN

      irrep_free

SOURCE

265 subroutine get_classes(nsym,sym,nclass,nelements,elements_idx)
266 
267 
268 !This section has been created automatically by the script Abilint (TD).
269 !Do not modify the following lines by hand.
270 #undef ABI_FUNC
271 #define ABI_FUNC 'get_classes'
272 !End of the abilint section
273 
274  implicit none
275 
276 !Arguments ------------------------------------
277 !scalars
278  integer,intent(in) :: nsym
279  integer,intent(out) :: nclass
280 !arrays
281  integer,intent(in) :: sym(3,3,nsym)
282  integer,intent(out) :: nelements(nsym),elements_idx(nsym,nsym)
283 
284 !Local variables-------------------------------
285 !scalars
286  integer :: isym,jsym,ksym,identity_idx !,ierr
287  character(len=500) :: msg
288 !arrays
289  integer :: cjg(3,3),ss(3,3),xx(3,3),xxm1(3,3),test(3,3)
290  integer :: identity(3,3)
291  logical :: found(nsym),found_identity
292 
293 !************************************************************************
294 
295  ! === Check if identity is present in the first position ===
296  identity=RESHAPE((/1,0,0,0,1,0,0,0,1/),(/3,3/)); found_identity=.FALSE.
297 
298  do isym=1,nsym
299    if (ALL(sym(:,:,isym)==identity)) then
300      found_identity=.TRUE.; identity_idx=isym; EXIT
301    end if
302  end do
303  if (.not.found_identity.or.identity_idx/=1) then
304   write(msg,'(3a)')&
305 &  'Either identity is not present or it is not the first operation ',ch10,&
306 &  'check set of symmetry operations '
307   MSG_ERROR(msg)
308  end if
309  !
310  ! Is it a group? Note that I assume that AFM sym.op (if any) have been pruned in the caller.
311  !dummy_symafm=1
312  !call chkgrp(nsym,dummy_symafm,sym,ierr)
313  !ABI_CHECK(ierr==0,"Error in group closure")
314 
315  nclass=0; nelements(:)=0; elements_idx(:,:)=0; found(:)=.FALSE.
316  do isym=1,nsym
317    if (.not.found(isym)) then
318      nclass=nclass+1
319      ss(:,:)=sym(:,:,isym)
320 
321      do jsym=1,nsym ! Form conjugate.
322        xx(:,:)=sym(:,:,jsym)
323        call mati3inv(xx,xxm1) ; xxm1=TRANSPOSE(xxm1)
324        cjg(:,:)=MATMUL(xxm1,MATMUL(ss,xx))
325        do ksym=1,nsym ! Is it already found?
326          test(:,:)=sym(:,:,ksym)
327          if (.not.found(ksym).and.(ALL((test-cjg)==0))) then
328            found(ksym)=.TRUE.
329            nelements(nclass)=nelements(nclass)+1
330            elements_idx(nelements(nclass),nclass)=ksym
331          end if
332        end do
333      end do
334 
335    end if
336  end do
337 
338 end subroutine get_classes

m_ptgroups/get_point_group [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

  get_point_group

FUNCTION

INPUTS

 ptg_name=point group name as returned by symptgroup.

OUTPUT

 nsym=Number of symmetries in the point group.
 nclass=Number of classes.
 sym(3,3,nsym)=Elements of the point group ordered by classe.
 class_ids(2,nclass)=Initial and final index in sym, for each
 Irreps(nclass)=Datatype gathering data on the different irreducible representations.

PARENTS

      m_ptgroups

CHILDREN

      irrep_free

SOURCE

121 subroutine get_point_group(ptg_name,nsym,nclass,sym,class_ids,class_names,Irreps)
122 
123 
124 !This section has been created automatically by the script Abilint (TD).
125 !Do not modify the following lines by hand.
126 #undef ABI_FUNC
127 #define ABI_FUNC 'get_point_group'
128 !End of the abilint section
129 
130  implicit none
131 
132 !Arguments ------------------------------------
133 !scalars
134  integer,intent(out) :: nclass,nsym
135  character(len=5),intent(in) :: ptg_name
136 !arrays
137  integer,allocatable,intent(out) :: sym(:,:,:),class_ids(:,:)
138  character(len=5),allocatable,intent(out) :: class_names(:)
139  type(irrep_t),allocatable,intent(out) :: Irreps(:)
140 
141 !Local variables-------------------------------
142 !scalars
143  integer :: irp,isym
144  !character(len=500) :: msg
145 
146 ! *************************************************************************
147 
148  SELECT CASE (TRIM(ADJUSTL(ptg_name)))
149  CASE ('1')
150    call ptg_C1  (nsym,nclass,sym,class_ids,class_names,Irreps)
151  CASE ('-1')
152    call ptg_Ci  (nsym,nclass,sym,class_ids,class_names,Irreps)
153  CASE ('2')
154    call ptg_C2  (nsym,nclass,sym,class_ids,class_names,Irreps)
155  CASE ('m',"-2") ! Abinit uses "-2"
156    call ptg_Cs  (nsym,nclass,sym,class_ids,class_names,Irreps)
157  CASE ('2/m')
158    call ptg_C2h (nsym,nclass,sym,class_ids,class_names,Irreps)
159  CASE ('222')
160    call ptg_D2  (nsym,nclass,sym,class_ids,class_names,Irreps)
161  CASE ('mm2')
162    call ptg_C2v (nsym,nclass,sym,class_ids,class_names,Irreps)
163  CASE ('mmm')
164    call ptg_D2h (nsym,nclass,sym,class_ids,class_names,Irreps)
165  CASE ('4')
166    call ptg_C4  (nsym,nclass,sym,class_ids,class_names,Irreps)
167  CASE ('-4')
168    call ptg_S4  (nsym,nclass,sym,class_ids,class_names,Irreps)
169  CASE ('4/m')
170    call ptg_C4h (nsym,nclass,sym,class_ids,class_names,Irreps)
171  CASE ('422')
172    call ptg_D4  (nsym,nclass,sym,class_ids,class_names,Irreps)
173  CASE ('4mm')
174    call ptg_C4v (nsym,nclass,sym,class_ids,class_names,Irreps)
175  CASE ('-42m')
176    call ptg_D2d (nsym,nclass,sym,class_ids,class_names,Irreps)
177  CASE ('4/mmm')
178    call ptg_D4h (nsym,nclass,sym,class_ids,class_names,Irreps)
179  CASE ('3')
180    call ptg_C3  (nsym,nclass,sym,class_ids,class_names,Irreps)
181  CASE ('-3')
182    call ptg_C3i (nsym,nclass,sym,class_ids,class_names,Irreps)
183  CASE ('32')
184    call ptg_D3  (nsym,nclass,sym,class_ids,class_names,Irreps)
185  CASE ('3m')
186    call ptg_C3v (nsym,nclass,sym,class_ids,class_names,Irreps)
187  CASE ('-3m')
188    call ptg_D3d (nsym,nclass,sym,class_ids,class_names,Irreps)
189  CASE ('6')
190    call ptg_C6  (nsym,nclass,sym,class_ids,class_names,Irreps)
191  CASE ('-6')
192    call ptg_C3h (nsym,nclass,sym,class_ids,class_names,Irreps)
193  CASE ('6/m')
194    call ptg_C6h (nsym,nclass,sym,class_ids,class_names,Irreps)
195  CASE ('622')
196    call ptg_D6  (nsym,nclass,sym,class_ids,class_names,Irreps)
197  CASE ('6mm')
198    call ptg_C6v (nsym,nclass,sym,class_ids,class_names,Irreps)
199  CASE ('-62m')
200    call ptg_D3h (nsym,nclass,sym,class_ids,class_names,Irreps)
201  CASE ('6/mmm')
202    call ptg_D6h (nsym,nclass,sym,class_ids,class_names,Irreps)
203  CASE ('23')
204    call ptg_T   (nsym,nclass,sym,class_ids,class_names,Irreps)
205  CASE ('m-3')
206    call ptg_Th  (nsym,nclass,sym,class_ids,class_names,Irreps)
207  CASE ('432')
208    call ptg_O   (nsym,nclass,sym,class_ids,class_names,Irreps)
209  CASE ('-43m')
210    call ptg_Td  (nsym,nclass,sym,class_ids,class_names,Irreps)
211  CASE ('m-3m')
212    call ptg_Oh  (nsym,nclass,sym,class_ids,class_names,Irreps)
213  CASE DEFAULT
214    MSG_BUG(sjoin("Unknown value for ptg_name:", ptg_name))
215  END SELECT
216 
217  ! Calculate the trace of each irreducible representation in order to have the character at hand.
218  do irp=1,SIZE(Irreps)
219    ABI_MALLOC(Irreps(irp)%trace, (nsym))
220    do isym=1,nsym
221      Irreps(irp)%trace(isym) = get_trace(Irreps(irp)%mat(:,:,isym))
222    end do
223  end do
224 
225 end subroutine get_point_group

m_ptgroups/groupk_free [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

  groupk_free

FUNCTION

  Deallocate all memory allocate in the group_k_t.

PARENTS

CHILDREN

      irrep_free

SOURCE

1241 subroutine groupk_free(Gk)
1242 
1243 
1244 !This section has been created automatically by the script Abilint (TD).
1245 !Do not modify the following lines by hand.
1246 #undef ABI_FUNC
1247 #define ABI_FUNC 'groupk_free'
1248 !End of the abilint section
1249 
1250  implicit none
1251 
1252 !Arguments ------------------------------------
1253  type(group_k_t),intent(inout) :: Gk
1254 
1255 ! *************************************************************************
1256 
1257 ! integer
1258  if (allocated(Gk%class_ids))   then
1259    ABI_FREE(Gk%class_ids)
1260  end if
1261  if (allocated(Gk%sym))   then
1262    ABI_FREE(Gk%sym)
1263  end if
1264 
1265 !real
1266  if (allocated(Gk%tnons)) then
1267    ABI_FREE(Gk%tnons)
1268  end if
1269 
1270 !character
1271  if (allocated(Gk%class_names)) then
1272    ABI_FREE(Gk%class_names)
1273  end if
1274 
1275 !type
1276  if (allocated(Gk%Irreps)) then
1277    call irrep_free(Gk%Irreps)
1278    ABI_DT_FREE(Gk%Irreps)
1279  end if
1280 
1281 end subroutine groupk_free

m_ptgroups/groupk_from_file [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

  groupk_from_file

FUNCTION

  Initialize the group_k_t datatype from an external database retrieved from
  the Bilbao server via the ptg.py script.

INPUTS

  fname(len=*)=file name

OUTPUT

  ierr=Status error
  Lgrps<group_k_t>=The structure completely initialized.

TODO

   This is a stub. I still have to complete the fileformat for the Bilbao database.

PARENTS

CHILDREN

      irrep_free

SOURCE

780 subroutine groupk_from_file(Lgrps,spgroup,fname,nkpt,klist,ierr)
781 
782 
783 !This section has been created automatically by the script Abilint (TD).
784 !Do not modify the following lines by hand.
785 #undef ABI_FUNC
786 #define ABI_FUNC 'groupk_from_file'
787 !End of the abilint section
788 
789  implicit none
790 
791 !Arguments ------------------------------------
792 !scalars
793  integer,intent(in) :: spgroup
794  integer,intent(out) :: ierr,nkpt
795  character(len=fnlen),intent(in) :: fname
796 !arrays
797  type(group_k_t),target,allocatable :: Lgrps(:)
798  real(dp),pointer :: klist(:,:)
799 
800 !Local variables-------------------------------
801 !scalars
802  integer,parameter :: last_file_version=1
803  integer :: unt,fvers,ik,nsym_ltgk,isym,nirreps_k,irp,icls
804  integer :: irrep_idx,irrep_dim,sym_idx,ita_spgroup,nels,prev,now
805  character(len=IRREPNAME_LEN) :: irrep_name
806  character(len=1) :: basis
807  character(len=500) :: msg
808  type(group_k_t),pointer :: Gk
809 !arrays
810  integer,allocatable :: nelements(:),elements_idx(:,:)
811  real(dp) :: kpt(3)
812  character(len=10),allocatable :: kname(:)
813  type(irrep_t),pointer :: OneIrr
814 
815 ! *************************************************************************
816 
817  ierr=0
818  if (open_file(fname,msg,newunit=unt,form="formatted") /=0) then
819    MSG_ERROR(msg)
820  end if
821 
822  read(unt,*,ERR=10)          ! Skip the header.
823  read(unt,*,ERR=10) fvers    ! File version.
824  if (fvers > last_file_version) then
825    write(msg,"(2(a,i0))")" Found file format= ",fvers," but the latest supported version is: ",last_file_version
826    MSG_ERROR(msg)
827  end if
828  read(unt,*,ERR=10) ita_spgroup
829  read(unt,*,ERR=10) basis
830 
831  if (spgroup/=ita_spgroup) then
832    write(msg,'(a,2i0)')&
833 &   " Input space group does not match with the value reported on file: ",spgroup,ita_spgroup
834    MSG_ERROR(msg)
835  end if
836 
837  if (basis /= "b") then
838    msg=" Wrong value for basis: "//TRIM(basis)
839    MSG_ERROR(msg)
840  end if
841 
842  ! * Read the list of the k-points.
843  read(unt,*,ERR=10) nkpt
844 
845  ABI_DT_MALLOC(Lgrps,(nkpt))
846 
847  ABI_MALLOC(klist,(3,nkpt))
848  ABI_MALLOC(kname,(nkpt))
849  do ik=1,nkpt
850    read(unt,*,ERR=10) klist(:,ik), kname(ik)
851  end do
852 
853  ! * Read tables for each k-point
854  do ik=1,nkpt
855 
856    read(unt,*,ERR=10) kpt
857    read(unt,*,ERR=10) nsym_ltgk
858    Gk  => Lgrps(ik)
859 
860    Gk%spgroup = ita_spgroup
861    Gk%nsym    = nsym_ltgk
862    Gk%point   = kpt
863    ABI_MALLOC(Gk%sym,(3,3,nsym_ltgk))
864    ABI_MALLOC(Gk%tnons,(3,nsym_ltgk))
865 
866    do isym=1,nsym_ltgk ! Read symmetries of the little group.
867      read(unt,*,ERR=10) Gk%sym(:,:,isym)
868      read(unt,*,ERR=10) Gk%tnons(:,isym)
869    end do
870 
871    ABI_MALLOC(nelements,(nsym_ltgk))
872    ABI_MALLOC(elements_idx,(nsym_ltgk,nsym_ltgk))
873 
874    call get_classes(nsym_ltgk,Gk%sym,Gk%nclass,nelements,elements_idx)
875 
876    ! The operations reported on the file are supposed to be packed in classes
877    ! otherwise one should perform a rearrangement of the indices.
878    prev = 0
879    do icls=1,Gk%nclass
880      do isym=1,nelements(icls)
881        now = elements_idx(isym,icls)
882        if ( (now-prev) /= 1 ) then
883          write(msg,"(2(a,i0))")&
884 &          " Symmetries on file are not ordered in classes. icls= ",icls,", isym= ",isym
885          MSG_ERROR(msg)
886        else
887          prev = now
888        end if
889      end do
890    end do
891 
892    ABI_MALLOC(Gk%class_ids,(2,Gk%nclass))
893    do icls=1,Gk%nclass
894      nels = nelements(icls)
895      Gk%class_ids(1,icls) = elements_idx(1,   icls)
896      Gk%class_ids(2,icls) = elements_idx(nels,icls)
897    end do
898 
899    ABI_FREE(nelements)
900    ABI_FREE(elements_idx)
901 
902    ! Read the irreducible representations.
903    read(unt,*,ERR=10) nirreps_k
904    ABI_CHECK(Gk%nclass == nirreps_k,"Gk%nclass /= nirreps_k")
905 
906    !$$ allocate(Gk%class_names(Gk%nclass))
907    ABI_DT_MALLOC(Gk%Irreps,(nirreps_k))
908 
909    do irp=1,nirreps_k
910      OneIrr =>  Gk%Irreps(irp)
911      read(unt,*,ERR=10) irrep_idx, irrep_dim, irrep_name
912      call init_irrep(OneIrr,nsym_ltgk,irrep_dim,irrep_name)
913      do isym=1,nsym_ltgk
914        read(unt,*,ERR=10) sym_idx, OneIrr%mat(:,:,isym)
915        ABI_CHECK(sym_idx==irp,"sym_idx/=irp!")
916        ! Matrix elements on file are in the form (rho, theta) with theta given in degrees.
917        call cmplx_sphcart(OneIrr%mat(:,:,isym),from="Sphere",units="Degrees")
918        OneIrr%trace(isym) = get_trace(OneIrr%mat(:,:,isym))
919      end do
920    end do
921 
922  end do
923 
924  close(unt)
925  RETURN
926  !
927  ! Handle IO-error.
928 10 ierr=1
929  close(unt)
930  RETURN
931 
932 end subroutine groupk_from_file

m_ptgroups/init_irrep [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

  alloc_irrep

FUNCTION

  Initialize an instance of the irrep_t datatype.

INPUTS

  nsym=The number of symmetries.
  irr_dim=The dimension of the irrep.
  [irr_name]=The name of theirrep. "???" is used if not given

OUTPUT

  Irrep<irrep_t>=

PARENTS

CHILDREN

SOURCE

1116 subroutine init_irrep(Irrep,nsym,irr_dim,irr_name)
1117 
1118 
1119 !This section has been created automatically by the script Abilint (TD).
1120 !Do not modify the following lines by hand.
1121 #undef ABI_FUNC
1122 #define ABI_FUNC 'init_irrep'
1123 !End of the abilint section
1124 
1125  implicit none
1126 
1127 !Arguments ------------------------------------
1128 !scalars
1129  integer,intent(in) :: nsym
1130 !arrays
1131  integer,intent(in) :: irr_dim
1132  character(len=*),optional,intent(in) :: irr_name
1133  type(irrep_t),intent(inout) :: Irrep
1134 
1135 !Local variables-------------------------------
1136  !character(len=500) :: msg
1137 ! *********************************************************************
1138 
1139  !@irrep_t
1140  Irrep%dim      = irr_dim
1141  Irrep%nsym     = nsym
1142  Irrep%name     = "???"
1143  if (PRESENT(irr_name)) Irrep%name = irr_name
1144 
1145  ABI_MALLOC(Irrep%mat,(irr_dim,irr_dim,nsym))
1146  Irrep%mat = czero
1147 
1148  ABI_MALLOC(Irrep%trace,(nsym))
1149  Irrep%trace = czero
1150 
1151 end subroutine init_irrep

m_ptgroups/irrep_free_0d [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

 irrep_free_0d

FUNCTION

  Deallocate all memory allocated in the irrep_t datatype.

PARENTS

      m_ptgroups

CHILDREN

      irrep_free

SOURCE

952 subroutine irrep_free_0d(Irrep)
953 
954 
955 !This section has been created automatically by the script Abilint (TD).
956 !Do not modify the following lines by hand.
957 #undef ABI_FUNC
958 #define ABI_FUNC 'irrep_free_0d'
959 !End of the abilint section
960 
961  implicit none
962 
963 !Arguments ------------------------------------
964  type(irrep_t),intent(inout) :: Irrep
965 
966 ! *********************************************************************
967 
968  !@irrep_t
969  if (allocated(Irrep%trace))  then
970    ABI_FREE(Irrep%trace)
971  end if
972  if (allocated(Irrep%mat))  then
973    ABI_FREE(Irrep%mat)
974  end if
975 
976 end subroutine irrep_free_0d

m_ptgroups/irrep_free_1d [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

 irrep_free_1d

FUNCTION

  Deallocate all memory allocated in the irrep_t datatype.

PARENTS

CHILDREN

      irrep_free

SOURCE

 995 subroutine irrep_free_1d(Irrep)
 996 
 997 
 998 !This section has been created automatically by the script Abilint (TD).
 999 !Do not modify the following lines by hand.
1000 #undef ABI_FUNC
1001 #define ABI_FUNC 'irrep_free_1d'
1002 !End of the abilint section
1003 
1004  implicit none
1005 
1006 !Arguments ------------------------------------
1007  type(irrep_t),intent(inout) :: Irrep(:)
1008 
1009 !Local variables-------------------------------
1010  integer :: irp
1011 ! *********************************************************************
1012 
1013  do irp=1,SIZE(Irrep)
1014    call irrep_free_0d(Irrep(irp))
1015  end do
1016 
1017 end subroutine irrep_free_1d

m_ptgroups/locate_sym [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

  locate_sym

FUNCTION

  Given a symmetry operation asym, this routine returns its index in the Ptg%sym
  array and the index of the class it belongs to.

INPUTS

OUTPUT

PARENTS

      m_esymm

CHILDREN

      irrep_free

SOURCE

621 subroutine locate_sym(Ptg,asym,sym_idx,cls_idx,ierr)
622 
623 
624 !This section has been created automatically by the script Abilint (TD).
625 !Do not modify the following lines by hand.
626 #undef ABI_FUNC
627 #define ABI_FUNC 'locate_sym'
628 !End of the abilint section
629 
630  implicit none
631 
632 !Arguments ------------------------------------
633 !scalars
634  integer,intent(out) :: sym_idx,cls_idx
635  integer,optional,intent(out) :: ierr
636  type(point_group_t),intent(in) :: Ptg
637 !arrays
638  integer,intent(in) :: asym(3,3)
639 
640 !Local variables-------------------------------
641 !scalars
642  integer :: isym,icls
643  character(len=500) :: msg
644 
645 ! *********************************************************************
646 
647  sym_idx = 0
648  do isym=1,Ptg%nsym
649    if (ALL(asym == Ptg%sym(:,:,isym) )) then
650      sym_idx = isym
651      EXIT
652    end if
653  end do
654 
655  cls_idx = 0
656  do icls=1,Ptg%nclass
657    if (sym_idx >= Ptg%class_ids(1,icls) .and. &
658 &      sym_idx <= Ptg%class_ids(2,icls) ) then
659      cls_idx = icls
660      EXIT
661    end if
662  end do
663 
664  if (PRESENT(ierr)) ierr=0
665  if (sym_idx==0 .or. cls_idx==0) then
666    write(msg,'(a,9(i0,1x),3a,i1,a,i1)')&
667 &    " Symmetry: ",asym," not found in point group table ",ch10,&
668 &    " sym_idx= ",sym_idx, " and cls_idx= ",cls_idx
669    if (PRESENT(ierr)) then
670      ierr=1
671      MSG_WARNING(msg)
672    else
673      MSG_ERROR(msg)
674    end if
675  end if
676 
677 end subroutine locate_sym

m_ptgroups/mult_table [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

 mult_table

FUNCTION

  Given a set of nsym 3x3 operations which are supposed to form a group,
  this routine constructs the multiplication table of the group.

INPUTS

 nsym=number of symmetry operation
 sym(3,3,nsym)=the operations

OUTPUT

  mtab(nsym,nsym)=The index of the product S_i * S_j in the input set sym.

PARENTS

CHILDREN

      irrep_free

SOURCE

704 subroutine mult_table(nsym,sym,mtab)
705 
706 
707 !This section has been created automatically by the script Abilint (TD).
708 !Do not modify the following lines by hand.
709 #undef ABI_FUNC
710 #define ABI_FUNC 'mult_table'
711 !End of the abilint section
712 
713  implicit none
714 
715 !Arguments ------------------------------------
716 !scalars
717  integer,intent(in) :: nsym
718 !arrays
719  integer,intent(in) :: sym(3,3,nsym)
720  integer,intent(out) :: mtab(nsym,nsym)
721 
722 !Local variables-------------------------------
723 !scalars
724  integer :: isym,jsym,ksym
725  !character(len=500) :: msg
726 !arrays
727  integer :: prod_ij(3,3),found(nsym)
728 
729 !************************************************************************
730 
731  do jsym=1,nsym
732    found(:)=0 ! Each symmetry should compare only once in a given (row|col).
733 
734    do isym=1,nsym
735      prod_ij = MATMUL(sym(:,:,isym),sym(:,:,jsym))
736      do ksym=1,nsym
737        if ( ALL(prod_ij == sym(:,:,ksym)) ) then
738          found(ksym)=found(ksym)+1
739          mtab(isym,jsym) = ksym
740        end if
741      end do
742    end do ! jsym
743 
744    if ( ANY(found /= 1)) then
745      write(std_out,*)"found = ",found
746      MSG_ERROR("Input elements do not form a group")
747    end if
748  end do ! isym
749 
750 end subroutine mult_table

m_ptgroups/point_group_free [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

 point_group_free

FUNCTION

  Deallocate all memory allocated in the point_group_t datatype.

PARENTS

      m_ptgroups

CHILDREN

      irrep_free

SOURCE

419 subroutine point_group_free(Ptg)
420 
421 
422 !This section has been created automatically by the script Abilint (TD).
423 !Do not modify the following lines by hand.
424 #undef ABI_FUNC
425 #define ABI_FUNC 'point_group_free'
426 !End of the abilint section
427 
428  implicit none
429 
430 !Arguments ------------------------------------
431  type(point_group_t),intent(inout) :: Ptg
432 
433 !Local variables-------------------------------
434 ! *********************************************************************
435 
436  !@point_group_t
437  if (allocated(Ptg%class_ids))   then
438    ABI_FREE(Ptg%class_ids)
439  end if
440  if (allocated(Ptg%sym)) then
441    ABI_FREE(Ptg%sym)
442  end if
443  if (allocated(Ptg%class_names)) then
444    ABI_FREE(Ptg%class_names)
445  end if
446 
447  if (allocated(Ptg%Irreps)) then
448    call irrep_free(Ptg%Irreps)
449    ABI_DT_FREE(Ptg%Irreps)
450  end if
451 
452 end subroutine point_group_free

m_ptgroups/point_group_init [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

 point_group_init

FUNCTION

  Creation method for the point_group_t datatype.

INPUTS

  ptg_name=The name of the point group (International conventions).

OUTPUT

  The datatype completely initialized.

PARENTS

      m_esymm,m_ptgroups

CHILDREN

      irrep_free

SOURCE

479 subroutine point_group_init(Ptg,ptg_name)
480 
481 
482 !This section has been created automatically by the script Abilint (TD).
483 !Do not modify the following lines by hand.
484 #undef ABI_FUNC
485 #define ABI_FUNC 'point_group_init'
486 !End of the abilint section
487 
488  implicit none
489 
490 !Arguments ------------------------------------
491 !scalars
492  character(len=5),intent(in) :: ptg_name
493  type(point_group_t),intent(inout) :: Ptg
494 ! *********************************************************************
495 
496  !@point_group_t
497  !call wrtout(std_out," Retrieving point group data for: "//TRIM(ptg_name),"COLL")
498 
499  Ptg%gname = ptg_name
500  call get_point_group(Ptg%gname,Ptg%nsym,Ptg%nclass,Ptg%sym,Ptg%class_ids,Ptg%class_names,Ptg%Irreps)
501 
502 end subroutine point_group_init

m_ptgroups/point_group_print [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

 point_group_print

FUNCTION

INPUTS

OUTPUT

PARENTS

      m_ptgroups

CHILDREN

      irrep_free

SOURCE

525 subroutine point_group_print(Ptg,header,unit,mode_paral,prtvol)
526 
527 
528 !This section has been created automatically by the script Abilint (TD).
529 !Do not modify the following lines by hand.
530 #undef ABI_FUNC
531 #define ABI_FUNC 'point_group_print'
532 !End of the abilint section
533 
534  implicit none
535 
536 !Arguments ------------------------------------
537 !scalars
538  integer,optional,intent(in) :: unit,prtvol
539  character(len=4),optional,intent(in) :: mode_paral
540  character(len=*),optional,intent(in) :: header
541  type(point_group_t),target,intent(in) :: Ptg
542 
543 !Local variables-------------------------------
544  integer :: my_unt,my_prtvol,irp,icls,sidx
545  complex(dpc) :: trace
546  character(len=4) :: my_mode
547  character(len=500) :: msg
548  type(irrep_t),pointer :: Row
549 ! *********************************************************************
550 
551  my_unt   =std_out; if (PRESENT(unit      )) my_unt   =unit
552  my_prtvol=0      ; if (PRESENT(prtvol    )) my_prtvol=prtvol
553  my_mode  ='COLL' ; if (PRESENT(mode_paral)) my_mode  =mode_paral
554 
555  msg=' ==== Point Group Table ==== '
556  if (PRESENT(header)) msg=' ==== '//TRIM(ADJUSTL(header))//' ==== '
557  call wrtout(my_unt,msg,my_mode)
558 
559  write(std_out,*)REPEAT("=",80)
560  write(std_out,*)" Point group : ",TRIM(Ptg%gname)," Number of symmetries ",Ptg%nsym," Number of classes    ",Ptg%nclass
561 
562  write(std_out,"(a6)",advance="no")"Class "
563  do icls=1,Ptg%nclass
564    write(std_out,"('|',a10)",advance="no")Ptg%class_names(icls)
565  end do
566  write(std_out,"('|')",advance="no")
567  write(std_out,*)" "
568 
569  write(std_out,"(a6)",advance="no")"Mult  "
570  do icls=1,Ptg%nclass
571    write(std_out,"('|',i10)",advance="no")Ptg%class_ids(2,icls)-Ptg%class_ids(1,icls) + 1
572  end do
573  write(std_out,"('|')",advance="no")
574  write(std_out,*)" "
575 
576  do irp=1,SIZE(Ptg%Irreps)
577    Row =>  Ptg%Irreps(irp)
578    write(std_out,'(a6)',advance="no")TRIM(Row%name)
579 
580    do icls=1,Ptg%nclass
581      sidx = Ptg%class_ids(1,icls)
582      trace = Row%trace(sidx)
583      if (ABS(AIMAG(trace)) > 10**(-6)) then
584         write(std_out,"('|',(2f5.2))",advance="no")trace
585       else
586         write(std_out,"('|',(f10.2))",advance="no")REAL(trace)
587       end if
588    end do
589 
590    write(std_out,"('|')",advance="no")
591    write(std_out,*)" "
592  end do
593 
594  write(std_out,*)REPEAT("=",80)
595 
596 end subroutine point_group_print

m_ptgroups/show_character_tables [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

  show_character_tables

FUNCTION

   Printout of the caracter tables of the 32 point groups.

INPUTS

  [unit]=Unit number of output file. Defaults to std_out

OUTPUT

  Only writing.

PARENTS

CHILDREN

      irrep_free

SOURCE

363 subroutine show_character_tables(unit)
364 
365 
366 !This section has been created automatically by the script Abilint (TD).
367 !Do not modify the following lines by hand.
368 #undef ABI_FUNC
369 #define ABI_FUNC 'show_character_tables'
370 !End of the abilint section
371 
372  implicit none
373 
374 !Arguments ------------------------------------
375 !scalars
376  integer,optional,intent(in) :: unit
377 
378 !Local variables-------------------------------
379  integer :: igrp,my_unt
380  character(len=5) :: ptg_name
381  type(point_group_t) :: Ptg
382 !arrays
383  !integer,allocatable :: elements_idx(:,:),nelements(:)
384 
385 ! *********************************************************************
386 
387  my_unt = std_out; if (PRESENT(unit)) my_unt=unit
388 
389  do igrp=1,SIZE(ptgroup_names)
390    ptg_name = ptgroup_names(igrp)
391    call point_group_init(Ptg,ptg_name)
392    call point_group_print(Ptg,unit=my_unt)
393    !allocate(nelements(Ptg%nsym),elements_idx(Ptg%nsym,Ptg%nsym))
394    !call get_classes(Ptg%nsym,Ptg%sym,nclass,nelements,elements_idx)
395    !deallocate(nelements,elements_idx)
396    call point_group_free(Ptg)
397  end do
398 
399 end subroutine show_character_tables

m_ptgroups/sum_irreps [ Functions ]

[ Top ] [ m_ptgroups ] [ Functions ]

NAME

  sum_irreps

FUNCTION

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

1172 function sum_irreps(Irrep1,Irrep2,ii,jj,kk,ll) result(res)
1173 
1174 
1175 !This section has been created automatically by the script Abilint (TD).
1176 !Do not modify the following lines by hand.
1177 #undef ABI_FUNC
1178 #define ABI_FUNC 'sum_irreps'
1179 !End of the abilint section
1180 
1181  implicit none
1182 
1183 !Arguments ------------------------------------
1184 !scalars
1185  integer,intent(in) :: ii,jj,kk,ll
1186 !arrays
1187  type(irrep_t),intent(in) :: Irrep1,Irrep2
1188  complex(dpc) :: res
1189 
1190 !Local variables-------------------------------
1191  integer :: isym,nsym,ierr
1192  !character(len=500) :: msg
1193 ! *********************************************************************
1194 
1195  ierr=0; res = czero
1196 
1197  nsym = Irrep1%nsym
1198  if (nsym /= Irrep2%nsym) then
1199    MSG_WARNING("Irreps have different nsym")
1200    ierr=ierr+1
1201  end if
1202 
1203  if (Irrep1%dim /= Irrep2%dim) then
1204    MSG_WARNING("Irreps have different dimensions")
1205    write(std_out,*)Irrep1%dim,Irrep2%dim
1206    ierr=ierr+1
1207  end if
1208 
1209  if (ii>Irrep2%dim .or. jj>Irrep2%dim .or. &
1210 &    kk>Irrep1%dim .or. ll>Irrep1%dim) then
1211    MSG_WARNING("Wrong indeces")
1212    write(std_out,*)ii,Irrep2%dim,jj,Irrep2%dim,kk>Irrep1%dim,ll,Irrep1%dim
1213    ierr=ierr+1
1214  end if
1215 
1216  if (ierr/=0) RETURN
1217 
1218  do isym=1,nsym
1219    res = res + DCONJG(Irrep1%mat(ii,jj,isym)) * Irrep2%mat(kk,ll,isym)
1220  end do
1221 
1222 end function sum_irreps