TABLE OF CONTENTS


ABINIT/m_esymm [ Modules ]

[ Top ] [ Modules ]

NAME

 m_esymm

FUNCTION

 This module defines structures and provides procedures used to find
 the irreducible representations associated to electronic eigenstates.

COPYRIGHT

  Copyright (C) 2008-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 .

PARENTS

CHILDREN

SOURCE

22 #if defined HAVE_CONFIG_H
23 #include "config.h"
24 #endif
25 
26 #include "abi_common.h"
27 
28 MODULE m_esymm
29 
30  use defs_basis
31  use m_abicore
32  use m_errors
33 
34  use m_io_tools,       only : file_exists
35  use m_symtk,          only : matr3inv, chkgrp, symrelrot, littlegroup_q
36  use m_symfind,        only : symbrav
37  use m_fstrings,       only : int2char10, itoa, sjoin
38  use m_numeric_tools,  only : print_arr, set2unit, get_trace
39  use m_hide_lapack,    only : xgeev, xginv
40  use m_crystal,        only : crystal_t, idx_spatial_inversion
41  use m_defs_ptgroups,  only : point_group_t, irrep_t
42  use m_ptgroups,       only : get_classes, point_group_init, irrep_free,&
43 &                             copy_irrep, init_irrep, mult_table, sum_irreps
44 
45  implicit none
46 
47  private

m_esymm/esymm_failed [ Functions ]

[ Top ] [ m_esymm ] [ Functions ]

NAME

  esymm_failed

FUNCTION

INPUTS

  esymm<esymm_t>

PARENTS

SOURCE

1537 function esymm_failed(esymm)
1538 
1539 
1540 !This section has been created automatically by the script Abilint (TD).
1541 !Do not modify the following lines by hand.
1542 #undef ABI_FUNC
1543 #define ABI_FUNC 'esymm_failed'
1544 !End of the abilint section
1545 
1546  implicit none
1547 
1548 !Arguments ------------------------------------
1549 !scalars
1550  logical :: esymm_failed
1551  type(esymm_t),intent(in) :: esymm
1552 
1553 ! *********************************************************************
1554 
1555  esymm_failed = (esymm%err_status /= ESYM_NOERROR)
1556 
1557 end function esymm_failed

m_esymm/esymm_finalize [ Functions ]

[ Top ] [ m_esymm ] [ Functions ]

NAME

 esymm_finalize

FUNCTION

INPUTS

OUTPUT

PARENTS

      classify_bands

CHILDREN

      xgeev,xginv,zpotrf,ztrsm

SOURCE

1113 subroutine esymm_finalize(esymm,prtvol)
1114 
1115 
1116 !This section has been created automatically by the script Abilint (TD).
1117 !Do not modify the following lines by hand.
1118 #undef ABI_FUNC
1119 #define ABI_FUNC 'esymm_finalize'
1120 !End of the abilint section
1121 
1122  implicit none
1123 
1124 !Arguments ------------------------------------
1125 !scalars
1126  integer,intent(in) :: prtvol
1127  type(esymm_t),target,intent(inout) :: esymm
1128 
1129 !Local variables-------------------------------
1130  integer :: idg,ib1,ib2,idx,nunknown,dg_dim
1131  integer :: try,irep,nitems,nseen,isn
1132  integer :: isym,idg1,idg2,dim_mat,irr_idx2,irr_idx1
1133  real(dp),parameter :: TOL_TRACE=0.1_dp,TOL_ORTHO=0.1_dp,TOL_UNITARY=0.1_dp ! Large tolerance is needed to avoid problems.
1134  !real(dp),parameter :: TOL_TRACE=0.01_dp,TOL_ORTHO=0.01_dp,TOL_UNITARY=0.01_dp ! Large tolerance is needed to avoid problems.
1135  !real(dp),parameter :: TOL_TRACE=tol3,TOL_ORTHO=tol3,TOL_UNITARY=tol3 ! Large tolerance is needed to avoid problems.
1136  real(dp) :: uerr,max_err
1137  complex(dpc) :: ctest
1138  logical :: isnew
1139  character(len=500) :: msg
1140 !arrays
1141  integer,allocatable :: dims_seen(:)
1142  complex(dpc),allocatable :: traces_seen(:,:)
1143  complex(dpc),pointer :: trace(:)
1144  complex(dpc),pointer :: calc_mat(:,:),trace1(:),trace2(:)
1145  complex(dpc),allocatable :: cidentity(:,:)
1146 
1147 ! *************************************************************************
1148 
1149  !@esymm_t
1150 
1151  ! Each band is initialized as "Unknown".
1152  esymm%b2irrep = 0
1153 
1154  ! Force the matrices to be unitary.
1155  call polish_irreps(esymm%Calc_irreps)
1156 
1157  if (.not.esymm%has_chtabs) then
1158 
1159    write(msg,'(5a)')&
1160 &    "Reference character table not available. ",ch10,&
1161 &    "Symmetry analysis not available. Using heuristic method to classify the states.",ch10,&
1162 &    "It might not work, especially if accidental degeneracies are present."
1163    MSG_WARNING(msg)
1164    !
1165    ! The simplest thing we can do here is using the calculated matrices to get the
1166    ! character and comparing the results hoping everything is OK.
1167    ABI_MALLOC(traces_seen,(esymm%nsym_gk,esymm%ndegs))
1168    ABI_MALLOC(dims_seen,(esymm%ndegs))
1169 
1170    traces_seen=czero; nseen=1
1171    traces_seen(:,1) = esymm%Calc_irreps(1)%trace
1172    dims_seen(1)     = esymm%Calc_irreps(1)%dim
1173 
1174    do idg=2,esymm%ndegs
1175      dg_dim = esymm%Calc_irreps(idg)%dim
1176      trace => esymm%Calc_irreps(idg)%trace
1177      isnew=.TRUE.
1178      do isn=1,nseen
1179        if (ALL (ABS(trace - traces_seen(:,isn)) < TOL_TRACE) ) then
1180          isnew=.FALSE.; EXIT
1181        end if
1182      end do
1183 
1184      if (isnew) then
1185        nseen = nseen+1
1186        traces_seen(:,nseen) = trace
1187        dims_seen(nseen) = dg_dim
1188      end if
1189    end do
1190 
1191    if (nseen>esymm%nclass) then
1192      write(msg,'(3a)')&
1193 &      "The number of different calculated traces is found to be greater than nclasses!",ch10,&
1194 &      "Heuristic method clearly failed. Symmetry analysis cannot be performed."
1195      MSG_WARNING(msg)
1196      esymm%err_status = ESYM_HEUR_WRONG_NCLASSES
1197      esymm%err_msg    = msg
1198 
1199      do isn=1,nseen
1200        write(msg,'(a,i0)')" Representation: ",isn
1201        call wrtout(std_out,msg,"COLL")
1202        call print_arr(traces_seen(:,isn),max_r=esymm%nsym_gk,unit=std_out,mode_paral="COLL")
1203      end do
1204 
1205    else  ! It seems that the Heuristic method succeeded.
1206      do idg=1,esymm%ndegs
1207        ib1=esymm%degs_bounds(1,idg)
1208        ib2=esymm%degs_bounds(2,idg)
1209        trace => esymm%Calc_irreps(idg)%trace
1210        do isn=1,nseen
1211          if (ALL (ABS(trace - traces_seen(:,isn)) < TOL_TRACE) ) then
1212            esymm%b2irrep(ib1:ib2)=isn
1213            if (esymm%Calc_irreps(idg)%dim /= dims_seen(isn)) then
1214              write(msg,'(3a)')&
1215 &              "Found two set of degenerate states with same character but different dimension!",ch10,&
1216 &              "heuristic method clearly failed. Symmetry analysis cannot be performed."
1217              MSG_ERROR(msg)
1218              esymm%err_status = ESYM_HEUR_WRONG_DIMS
1219              esymm%err_msg    = msg
1220            end if
1221            EXIT
1222          end if
1223        end do
1224      end do
1225    end if
1226 
1227    ABI_FREE(traces_seen)
1228    ABI_FREE(dims_seen)
1229 
1230  else
1231    !
1232    ! * Search in the lookup table definining the irreducible representation
1233    nunknown = 0
1234    do idg=1,esymm%ndegs
1235 
1236      ib1=esymm%degs_bounds(1,idg)
1237      ib2=esymm%degs_bounds(2,idg)
1238      trace => esymm%Calc_irreps(idg)%trace
1239 
1240      try = which_irrep(esymm, trace, tol3)
1241      if (try==0) try = which_irrep(esymm, trace, 0.1_dp) ! try again with increased tolerance.
1242      if (try/=0) then
1243        esymm%b2irrep(ib1:ib2)=try
1244      else
1245        esymm%b2irrep(ib1:ib2)=0
1246        nunknown = nunknown + (ib2-ib1+1)
1247      end if
1248    end do
1249  end if
1250  !
1251  ! %irrep2b(0)) gives the indeces of the states that have not been classified.
1252  ABI_DT_MALLOC(esymm%irrep2b,(0:esymm%nclass))
1253 
1254  !write(std_out,*)"b2irrep",esymm%b2irrep
1255 
1256  do irep=0,esymm%nclass
1257    nitems = COUNT(esymm%b2irrep==irep)
1258    ABI_MALLOC(esymm%irrep2b(irep)%value,(nitems))
1259    idx=0
1260    do ib1=1,esymm%nbnds
1261      if (esymm%b2irrep(ib1) == irep) then
1262        idx = idx + 1
1263        esymm%irrep2b(irep)%value(idx) = ib1
1264      end if
1265    end do
1266  end do
1267 
1268  if (size(esymm%irrep2b(0)%value) /= 0) then
1269    write(msg,'(a,i0,a)')" Band classification algorithm was not able to classify ",size(esymm%irrep2b(0)%value)," states."
1270    MSG_WARNING(msg)
1271    esymm%err_status = ESYM_CLASSIFICATION_ERROR
1272    esymm%err_msg    = msg
1273  end if
1274  !
1275  ! ==============================================================
1276  ! ==== Test basic properties of irreducible representations ====
1277  ! ==============================================================
1278 
1279  if (.not.esymm_failed(esymm)) then
1280    !
1281    ! 1) \sum_R \chi^*_a(R)\chi_b(R)= N_R \delta_{ab}
1282    !
1283    !call wrtout(std_out," \sum_R \chi^*_a(R)\chi_b(R) = N_R \delta_{ab} ","COLL")
1284    max_err=zero
1285    do idg2=1,esymm%ndegs
1286      trace2 => esymm%Calc_irreps(idg2)%trace(1:esymm%nsym_gk)
1287      ib2 = esymm%degs_bounds(1,idg2)
1288      irr_idx2 = esymm%b2irrep(ib2)
1289      if (irr_idx2 == 0) CYCLE
1290 
1291      do idg1=1,idg2
1292        trace1 => esymm%Calc_irreps(idg1)%trace(1:esymm%nsym_gk)
1293        ib1 = esymm%degs_bounds(1,idg1)
1294        irr_idx1 = esymm%b2irrep(ib1)
1295        if (irr_idx1 == 0) CYCLE
1296        ctest=DOT_PRODUCT(trace1,trace2)/esymm%nsym_gk
1297        if (irr_idx1==irr_idx2) ctest=ctest-one
1298        max_err = MAX(max_err,ABS(ctest))
1299        if (.FALSE..and.ABS(ctest)>tol3) then
1300          write(msg,'(a,4i3,2es16.8)')&
1301 &          ' WARNING: should be delta_ij: cx1 cx2, irr1, irr2, ctest: ',idg1,idg2,irr_idx1,irr_idx2,ctest
1302          call wrtout(std_out,msg,"COLL")
1303        end if
1304      end do
1305    end do
1306 
1307    if (max_err>TOL_ORTHO) then
1308      write(msg,'(a,es10.2)')" Too large maximum error on \sum_R \chi^*_a(R)\chi_b(R) = N_R \delta_{ab}: ",max_err
1309      MSG_WARNING(msg)
1310      esymm%err_status =  ESYM_ORTHO_ERROR
1311      esymm%err_msg    =  msg
1312    else
1313      write(msg,'(a,es10.2)')" maximum error on \sum_R \chi^*_a(R)\chi_b(R) = N_R \delta_{ab}: ",max_err
1314      call wrtout(std_out,msg,"COLL")
1315    end if
1316 
1317    if (.not.esymm%only_trace) then
1318      !call wrtout(std_out," **** Testing the unitary of the calculated irreps ****",my_mode)
1319      max_err=zero
1320      do idg1=1,esymm%ndegs
1321        ib1 = esymm%degs_bounds(1,idg1)
1322        irr_idx1 = esymm%b2irrep(ib1)
1323        if (irr_idx1 == 0) CYCLE
1324 
1325        do isym=1,esymm%nsym_gk
1326          calc_mat => esymm%Calc_irreps(idg1)%mat(:,:,isym)
1327          dim_mat  =  esymm%Calc_irreps(idg1)%dim
1328          ABI_MALLOC(cidentity,(dim_mat,dim_mat))
1329          call set2unit(cidentity)
1330          uerr = MAXVAL( ABS(MATMUL(calc_mat,TRANSPOSE(DCONJG(calc_mat))) - cidentity) )
1331          max_err = MAX(max_err,uerr)
1332          ABI_FREE(cidentity)
1333          if (.FALSE..and.prtvol>=10) then
1334            write(std_out,'(a,i3,a,i2,a,es16.8,a)')&
1335 &          " === idg: ",idg1,", isym: ",isym,", Error on U^* U = 1: ",uerr," ==="
1336            call print_arr(calc_mat,dim_mat,dim_mat,unit=std_out,mode_paral="COLL")
1337          end if
1338        end do
1339      end do
1340 
1341      if (max_err>TOL_UNITARY) then
1342        write(msg,'(a,es10.2)')" Too large maximum error on the unitary of representions matrices: ",max_err
1343        MSG_WARNING(msg)
1344        esymm%err_msg    = msg
1345        esymm%err_status = ESYM_UNITARY_ERROR
1346      else
1347        write(msg,'(a,es10.2)')" maximum error on the unitary of representions matrices: ",max_err
1348        call wrtout(std_out,msg,"COLL")
1349      end if
1350 
1351    end if
1352 
1353  end if
1354 
1355 end subroutine esymm_finalize

m_esymm/esymm_free_0D [ Functions ]

[ Top ] [ m_esymm ] [ Functions ]

NAME

 esymm_free_0D

FUNCTION

  Deallocate the memory allocated in the esymm_t datatype (scalar version)

PARENTS

      m_esymm

CHILDREN

      xgeev,xginv,zpotrf,ztrsm

SOURCE

 978 subroutine esymm_free_0D(esymm)
 979 
 980 
 981 !This section has been created automatically by the script Abilint (TD).
 982 !Do not modify the following lines by hand.
 983 #undef ABI_FUNC
 984 #define ABI_FUNC 'esymm_free_0D'
 985 !End of the abilint section
 986 
 987  implicit none
 988 
 989 !Arguments ------------------------------------
 990 !scalars
 991  type(esymm_t),intent(inout) :: esymm
 992 
 993 !Local variables ------------------------------
 994 !scalars
 995  integer :: ii
 996 ! *************************************************************************
 997 
 998  !@esymm_t
 999  if (allocated(esymm%g0)) then
1000    ABI_FREE(esymm%g0)
1001  end if
1002  if (allocated(esymm%tr_g0))  then
1003    ABI_FREE(esymm%tr_g0)
1004  end if
1005  if (allocated(esymm%nelements))  then
1006    ABI_FREE(esymm%nelements)
1007  end if
1008  if (allocated(esymm%sgk2symrec))  then
1009    ABI_FREE(esymm%sgk2symrec)
1010  end if
1011  if (allocated(esymm%tr_sgk2symrec)) then
1012    ABI_FREE(esymm%tr_sgk2symrec)
1013  end if
1014  if (allocated(esymm%herring_test)) then
1015    ABI_FREE(esymm%herring_test)
1016  end if
1017  if (allocated(esymm%b2irrep)) then
1018    ABI_FREE(esymm%b2irrep)
1019  end if
1020  if (allocated(esymm%degs_bounds))  then
1021    ABI_FREE(esymm%degs_bounds)
1022  end if
1023  if (allocated(esymm%degs_dim)) then
1024    ABI_FREE(esymm%degs_dim)
1025  end if
1026 
1027  if (allocated(esymm%irrep2b)) then
1028    do ii=LBOUND(esymm%irrep2b,DIM=1),UBOUND(esymm%irrep2b,DIM=1)
1029      ABI_FREE(esymm%irrep2b(ii)%value)
1030    end do
1031    ABI_DT_FREE(esymm%irrep2b)
1032  end if
1033 
1034  if (allocated(esymm%Calc_irreps)) then
1035    call irrep_free(esymm%Calc_irreps)
1036  end if
1037 
1038  if (allocated(esymm%trCalc_irreps)) then
1039    call irrep_free(esymm%trCalc_irreps)
1040  end if
1041 
1042  if (allocated(esymm%Ref_irreps)) then
1043    call irrep_free(esymm%Ref_irreps)
1044  end if
1045 
1046 end subroutine esymm_free_0D

m_esymm/esymm_free_2D [ Functions ]

[ Top ] [ m_esymm ] [ Functions ]

NAME

 esymm_free_2D

FUNCTION

  Deallocate the memory allocated in the esymm_t datatype (2D version)

PARENTS

CHILDREN

      xgeev,xginv,zpotrf,ztrsm

SOURCE

1065 subroutine esymm_free_2D(esymm)
1066 
1067 
1068 !This section has been created automatically by the script Abilint (TD).
1069 !Do not modify the following lines by hand.
1070 #undef ABI_FUNC
1071 #define ABI_FUNC 'esymm_free_2D'
1072 !End of the abilint section
1073 
1074  implicit none
1075 
1076 !Arguments ------------------------------------
1077 !scalars
1078  type(esymm_t),intent(inout) :: esymm(:,:)
1079 
1080 !Local variables ------------------------------
1081  integer :: id1,id2
1082 ! *************************************************************************
1083 
1084  do id2=1,SIZE(esymm,DIM=2)
1085    do id1=1,SIZE(esymm,DIM=1)
1086      call esymm_free_0D(esymm(id1,id2))
1087    end do
1088  end do
1089 
1090 end subroutine esymm_free_2D

m_esymm/esymm_init [ Functions ]

[ Top ] [ m_esymm ] [ Functions ]

NAME

 esymm_init

FUNCTION

  Initialize a esymm_t datatype containing data and parameters
  needed to analyze the irreducible representations at a particular k-point....

INPUTS

  kpt_in(3)=The k-point where the classification of bands is required.
  Cryst<crystal_t>=Datatype describing the unit cell and its symmetries.
  nspinor=number of spinorial components
  nsppol=number of independent polarizations
  first_ib=Index of the first band.
  nbnds=Number of bands for this k-point.
  ene_k(nbnds)=energies for this k-point. ene_k(1) corresponds to band first_ib.
  EDIFF_TOL=tolerance below which two states are considered to belong to the same irreducible representation

OUTPUT

  esymm<esymm_t>= Initialized data type gathering information of the small group
     of the k-point as well as the irreducible representations.

NOTES

   The present implementation does NOT work at zone border if the little group of
   kpt_in is non-symmorphic namely thers is at lest a symmetry operation with non-zero tnons.

PARENTS

      classify_bands

CHILDREN

      xgeev,xginv,zpotrf,ztrsm

SOURCE

257 subroutine esymm_init(esymm,kpt_in,Cryst,only_trace,nspinor,first_ib,nbnds,EDIFF_TOL,ene_k,tolsym)
258 
259 
260 !This section has been created automatically by the script Abilint (TD).
261 !Do not modify the following lines by hand.
262 #undef ABI_FUNC
263 #define ABI_FUNC 'esymm_init'
264 !End of the abilint section
265 
266  implicit none
267 
268 !Arguments ------------------------------------
269 !scalars
270  integer,intent(in) :: nbnds,nspinor,first_ib
271  real(dp),intent(in) :: EDIFF_TOL,tolsym
272  logical,intent(in) :: only_trace
273  type(crystal_t),intent(in) :: Cryst
274  type(esymm_t),intent(out) :: esymm
275 !arrays
276  real(dp),intent(in) :: ene_k(nbnds),kpt_in(3)
277 
278 !Local variables-------------------------------
279 !scalars
280  integer :: dim_degs,iband,idg,irp,nacc_deg,isym_gk,grp_ierr
281  integer :: nsym_fm,idx_fm,idx_gk,idx_trgk,isym,jsym,dummy_timrev !,iholohedry
282  integer :: iel,icls,msym,iord !isym1,!iprod,dim_irrep,icls2, isym2,isym_tr,
283  integer :: spgroup,chkprim !,ptgroupma
284  real(dp) :: mkt
285  !complex(dpc) :: phase_k
286  character(len=5) :: ptgroup,ptgroup_name
287  character(len=10) :: spgroup_str
288  character(len=1000) :: msg
289  character(len=fnlen) :: lgroup_fname
290 !arrays
291  integer :: inversion(3,3)
292  integer,allocatable :: degs_bounds(:,:),dim_irreps(:)
293  integer :: bravais(11),sym_axis(3)
294  real(dp) :: pmat1(3,3),pmat2(3,3),pmat3(3,3),pmat4(3,3),pmat5(3,3),pmat6(3,3)
295  !real(dp) :: genafm(3)
296  !integer :: rot2(3,3)
297  !integer,allocatable :: mtab(:,:)
298  integer,allocatable :: elements_idx(:,:),tmp_nelements(:)
299  integer,allocatable :: found(:),symrec_fm(:,:,:),fm2symrec(:)
300  integer,allocatable :: ksym_table(:,:,:),sgk(:,:,:),tr_sgk(:,:,:),dum_symafm(:)
301  integer,allocatable :: new_idx(:),new_g0(:,:),tmp_symrec(:,:,:),conv_symrec(:,:,:) !,tr_conv_symrec(:,:,:)
302  integer,allocatable :: dummy_symafm(:)
303  real(dp) :: conv_gprimd(3,3),axes(3,3) !,tau2(3)
304  !complex(dpc),allocatable :: her_test(:) !,mat_test(:,:)
305  complex(dpc),allocatable :: phase_mkt(:)
306  type(point_group_t) :: Ptg
307 
308 ! *************************************************************************
309 
310  DBG_ENTER("COLL")
311 
312  !@esymm_t
313  esymm%err_status= ESYM_NOERROR
314  inversion=RESHAPE((/-1,0,0,0,-1,0,0,0,-1/),(/3,3/))
315  !
316  ! ====================================
317  ! ==== Initialize basic variables ====
318  ! ====================================
319  esymm%nspinor        = nspinor
320  esymm%first_ib       = first_ib
321  esymm%nbnds          = nbnds
322  esymm%only_trace     = only_trace
323  esymm%tol_deg        = EDIFF_TOL
324  esymm%has_spatial_inv= (idx_spatial_inversion(Cryst) /= 0)
325  esymm%can_use_tr     = .TRUE. !TODO this should be input
326  esymm%has_chtabs     = .FALSE.
327  esymm%kpt            = kpt_in(:)
328  esymm%nonsymmorphic_at_zoneborder=.FALSE.
329  !
330  ! ===============================
331  ! === Locate degenerate_bands ===
332  ! ===============================
333  esymm%ndegs=1
334 
335  ABI_MALLOC(degs_bounds,(2,nbnds))
336  degs_bounds=0; degs_bounds(1,1)=1
337 
338  do iband=2,nbnds
339    if (ABS(ene_k(iband)-ene_k(iband-1))>EDIFF_TOL) then
340      degs_bounds(2,esymm%ndegs) = iband-1 + (first_ib-1)
341      esymm%ndegs=esymm%ndegs+1
342      degs_bounds(1,esymm%ndegs) = iband + (first_ib-1)
343    end if
344  end do
345  degs_bounds(2,esymm%ndegs)=nbnds + (first_ib-1)
346 
347  ABI_MALLOC(esymm%degs_bounds,(2,esymm%ndegs))
348  esymm%degs_bounds = degs_bounds(:,1:esymm%ndegs)
349  ABI_FREE(degs_bounds)
350  !
351  ! Each band is initialized as "Unknown".
352  ABI_MALLOC(esymm%b2irrep,(esymm%nbnds))
353  esymm%b2irrep = 0
354  !
355  ! ==================================
356  ! ==== Find the group of kpt_in ====
357  ! ==================================
358  ! The small point group is the subset of symrec such that $ S q = q + g0 $
359  ! Symmetries are packed in classes.
360  ! For the time being, AFM symmetries are not treated.
361 
362  write(msg,'(a,3(1x,f7.4))')" Finding the little group of k-point: ",esymm%kpt
363  call wrtout(std_out,msg,"COLL")
364 
365  ! Only FM symmetries are used.
366  nsym_fm = COUNT(Cryst%symafm==1)
367 
368  if (nsym_fm /= Cryst%nsym) then
369    write(msg,'(4a)')ch10,&
370 &    "Band classification in terms of magnetic space groups not coded! ",ch10,&
371 &    "Only the ferromagnetic subgroup will be used "
372    MSG_COMMENT(msg)
373  end if
374 
375  ABI_MALLOC(symrec_fm,(3,3,nsym_fm))
376  ABI_MALLOC(fm2symrec,(nsym_fm))
377 
378  idx_fm = 0
379  do isym=1,Cryst%nsym
380    if (Cryst%symafm(isym) == 1) then
381      idx_fm = idx_fm+1
382      symrec_fm(:,:,idx_fm) = Cryst%symrec(:,:,isym)
383      fm2symrec(idx_fm) = isym
384    end if
385  end do
386 
387  ! Find symmetries that preserve k.
388  ABI_MALLOC(ksym_table,(4,2,nsym_fm))
389  ABI_MALLOC(dummy_symafm,(nsym_fm))
390 
391  dummy_symafm = 1
392  call littlegroup_q(nsym_fm,esymm%kpt,ksym_table,symrec_fm,dummy_symafm,dummy_timrev,prtvol=0)
393 
394  esymm%nsym_gk  =COUNT(ksym_table(4,1,:)==1)  ! # S such that  S k = k +G0
395 
396  esymm%nsym_trgk=0
397  if (esymm%can_use_tr) esymm%nsym_trgk=COUNT(ksym_table(4,2,:)==1)  ! # S such that -S k = k +G0
398 
399  ! Allocate workspace arrays.
400  ABI_MALLOC(sgk,(3,3,esymm%nsym_gk))
401  ABI_MALLOC(tr_sgk,(3,3,esymm%nsym_trgk))
402 
403  ! Allocate mapping little-group --> symrec and table for umklapps.
404  ABI_MALLOC(esymm%sgk2symrec,(esymm%nsym_gk))
405  ABI_MALLOC(esymm%g0,(3,esymm%nsym_gk))
406  ABI_MALLOC(esymm%tr_sgk2symrec,(esymm%nsym_trgk))
407  ABI_MALLOC(esymm%tr_g0,(3,esymm%nsym_trgk))
408 
409  ! Important NOTE:
410  ! If nonsymmorphic_at_zoneborder symmetry analysis cannot be performed unless
411  ! an external database retrieved from the bilbao server (REPRES) is found.
412  idx_gk=0; idx_trgk=0
413  esymm%sgk2symrec=-999; esymm%tr_sgk2symrec=-999
414  do isym=1,nsym_fm
415 
416    if (ksym_table(4,1,isym)==1) then ! S k = k +G0
417      idx_gk=idx_gk+1
418      sgk(:,:,idx_gk)=symrec_fm(:,:,isym)
419      esymm%g0(:,idx_gk)=ksym_table(1:3,1,isym)
420      esymm%sgk2symrec(idx_gk)=fm2symrec(isym)
421      if (ANY(ksym_table(1:3,1,isym)/=0).and.(ANY(ABS(Cryst%tnons(:,fm2symrec(isym)))>tol6))) then
422         esymm%nonsymmorphic_at_zoneborder=.TRUE.
423      end if
424    end if
425 
426    if (esymm%can_use_tr.and.ksym_table(4,2,isym)==1) then ! -S k = k +G0
427      idx_trgk=idx_trgk+1
428      tr_sgk(:,:,idx_trgk)=symrec_fm(:,:,isym)
429      esymm%tr_g0(:,idx_trgk)=ksym_table(1:3,2,isym)
430      esymm%tr_sgk2symrec(idx_trgk)=fm2symrec(isym)
431    end if
432  end do
433 
434  ABI_FREE(ksym_table)
435  ABI_FREE(symrec_fm)
436  ABI_FREE(fm2symrec)
437 
438 ! ==========================================
439 ! ==== Divide the operations in classes ====
440 ! ==========================================
441  ABI_MALLOC(dum_symafm,(esymm%nsym_gk))
442  dum_symafm=1
443 
444  call chkgrp(esymm%nsym_gk,dum_symafm,sgk,grp_ierr)
445 
446  ABI_CHECK(grp_ierr==0,"chkgrp failed")
447  ABI_FREE(dum_symafm)
448 
449  ABI_MALLOC(tmp_nelements,(esymm%nsym_gk))
450  ABI_MALLOC(elements_idx,(esymm%nsym_gk,esymm%nsym_gk))
451 
452  call get_classes(esymm%nsym_gk,sgk,esymm%nclass,tmp_nelements,elements_idx)
453 
454  ABI_MALLOC(esymm%nelements,(esymm%nclass))
455  esymm%nelements = tmp_nelements(1:esymm%nclass)
456  ABI_FREE(tmp_nelements)
457 
458  ! From the list of symmetry operations and the lattice vectors, determine the
459  ! Bravais information including the holohedry, the centering, the coordinate of
460  ! the primitive vectors in the conventional vectors, as well as the point group,
461  msym=192; if (allocated(Cryst%symrec)) msym=size(Cryst%symrec,3)
462  ABI_MALLOC(tmp_symrec,(3,3,msym))
463  tmp_symrec(:,:,1:esymm%nsym_gk)=sgk
464 
465  call symbrav(bravais,msym,esymm%nsym_gk,ptgroup,Cryst%gprimd,tmp_symrec,tolsym,axis=sym_axis)
466 
467  ABI_FREE(tmp_symrec)
468 
469  write(std_out,'(a)')" symptgroup returned point group: "//TRIM(ptgroup)
470  write(std_out,'(a,i2)')" iholohedry ",bravais(1)
471  write(std_out,'(a,i2)')" center     ",bravais(2)
472  write(std_out,'(a,9i3)')" gprimd in the axes of the conventional bravais lattice (*2 if center/=0)",bravais(3:11)
473  write(std_out,'(a,3i3)')" sym_axis ",sym_axis
474 
475  ! Branching:
476  ! 1) If the little group is not symmorphic_at_zoneborder we can
477  !    classify the states using the irreducible representation of the point group.
478  !
479  ! 2) If the little group is symmorphic_at_zoneborder, we have to rely on
480  !    an external database retrieved from the Bilbao server in order to classify the states.
481  !    If the file is not available, we only know the number of classes but neither their
482  !    character nor the dimension of the irreducible representation.
483  !
484  if (esymm%nonsymmorphic_at_zoneborder) then
485 
486    spgroup=0
487    chkprim=1 ! Cell must be primitive.
488    !call symlatt(bravais,msym,nptsym,ptsymrel,rprimd,tolsym)
489    !call symspgr(bravais,Cryst%nsym,spgroup,Cryst%symrel,Cryst%tnons,tolsym)
490 
491    !call symanal(bravais,chkprim,genafm,msym,nsym,ptgroupma,rprimd,spgroup,symafm,symrel,tnons,tolsym)
492 
493    call int2char10(spgroup,spgroup_str)
494    lgroup_fname = "lgroup_"//TRIM(spgroup_str)
495 
496    if (file_exists(lgroup_fname)) then
497      MSG_ERROR("Not coded")
498 
499      ! Read little groups from the external database.
500      !% call init_groupk_from_file(Lgrp,spgroup,lgroup_fname,ierr)
501 
502      ! Save the irreducible representations in esymm.
503      ! Reorder symmetries such that they correspond to the Bilbao database.
504      !% allocate(esymm%Ref_irreps(esymm%nclass))
505      !% call copy_irrep(Irreps, esymm%Ref_irreps)
506 
507    else
508      write(msg,'(7a)')&
509 &      "Non-symmorphic small group and zone border. ",ch10,&
510 &      "External file: ",TRIM(lgroup_fname)," containing Bilbao tables not found ",ch10,&
511 &      "Character analysis cannot be performed. Accidental degeneracies cannot be detected. "
512      MSG_WARNING(msg)
513 
514      esymm%has_chtabs = .FALSE.
515 
516      ! Reorder indeces such that symmetries are packed in classes.
517      ABI_MALLOC(new_idx,(esymm%nsym_gk))
518      ABI_MALLOC(new_g0,(3,esymm%nsym_gk))
519      new_g0=0; iord = 0
520      do icls=1,esymm%nclass
521        do iel=1,esymm%nelements(icls)
522          iord = iord+1
523          jsym = elements_idx(iel,icls)
524          new_idx(iord)  = esymm%sgk2symrec(jsym)
525          new_g0(:,iord) = esymm%g0(:,jsym)
526        end do
527      end do
528 
529      esymm%sgk2symrec = new_idx
530      esymm%g0 = new_g0
531 
532      ABI_FREE(new_idx)
533      ABI_FREE(new_g0)
534    end if ! file exists
535 
536  else
537    !
538    ! **** This part is still under development. It might not work for particular ****
539    ! **** orientations of the unit cell or particular lattices.                  ****
540    !
541    ! The symmetries in the Bilbao database refer to the conventional unit cells.
542    ! Therefore we have to map the abinit symmetries (in reduced coordinates)
543    ! onto the Bilbao dataset. Bilbao standard settings are:
544    !
545    ! * unique axis b (cell choice 1) for space groups withing the monoclinic system
546    ! * obverse triple hexagonal unit cell R space groups.
547    ! * origin choice two - inversion center at (0, 0, 0) - for the centrosymmetric
548    !   space groups for which there are two origins choices, within the
549    !   orthorombic, tetragonal and cubic system.
550 
551    ! 1) Retrieve the rotation matrices and the irreducible representations (Bilbao setting).
552    call point_group_init(Ptg,ptgroup)
553 
554    esymm%has_chtabs = .TRUE.
555    ABI_CHECK(esymm%nclass==Ptg%nclass,"esymm%nclass/=Ptg%nclass!")
556 
557    do icls=1,esymm%nclass ! FIXME this is awful, should be done in a cleaner way.
558      esymm%nelements(icls)=Ptg%class_ids(2,icls) - Ptg%class_ids(1,icls) + 1
559    end do
560 
561    ! 2) Generate the symmetry operations in the conventional vector coordinates.
562    conv_gprimd(:,1)=bravais(3:5)
563    conv_gprimd(:,2)=bravais(6:8)
564    conv_gprimd(:,3)=bravais(9:11)
565 
566    axes = conv_gprimd
567    call matr3inv(conv_gprimd,axes) !; axes=TRANSPOSE(axes)
568 
569    conv_gprimd=MATMUL(Cryst%gprimd,TRANSPOSE(axes))
570    !conv_gprimd=MATMUL(axes,Cryst%gprimd)
571    !conv_gprimd=MATMUL(TRANSPOSE(axes),Cryst%gprimd)
572    !write(std_out,*)"conv_gprimd:", conv_gprimd
573 
574    ptgroup_name = ADJUSTL(ptgroup)
575 
576    select case (ptgroup_name)
577 
578    case ("3m","-3m")
579      call wrtout(std_out," Changing the conventional cell: rhombohedral --> triple hexagonal","COLL")
580      ! Transformation matrices: primitive rhombohedral --> triple hexagonal cell obverse setting. Table 5.1.3.1 ITA page 81.
581      pmat1 = RESHAPE( (/ 1,-1, 0, 0, 1,-1, 1, 1, 1/), (/3,3/) ) ! R1
582      pmat2 = RESHAPE( (/ 0, 1,-1,-1, 0, 1, 1, 1, 1/), (/3,3/) ) ! R2
583      pmat3 = RESHAPE( (/-1, 0, 1, 1,-1, 0, 1, 1, 1/), (/3,3/) ) ! R3
584      pmat4 = RESHAPE( (/-1, 1, 0, 0,-1, 1, 1, 1, 1/), (/3,3/) ) ! R1 reverse setting.
585      pmat5 = RESHAPE( (/ 0,-1, 1, 1, 0,-1, 1, 1, 1/), (/3,3/) ) ! R2 reverse setting.
586      pmat6 = RESHAPE( (/ 1, 0,-1,-1, 1, 0, 1, 1, 1/), (/3,3/) ) ! R3 reverse setting.
587      conv_gprimd = MATMUL(conv_gprimd,pmat1)
588      !conv_gprimd = MATMUL(conv_gprimd,pmat2)
589      !conv_gprimd = MATMUL(conv_gprimd,pmat3)
590      !conv_gprimd = MATMUL(conv_gprimd,pmat4)
591      !conv_gprimd = MATMUL(conv_gprimd,pmat5)
592      !conv_gprimd = MATMUL(conv_gprimd,pmat6)
593      !write(std_out,*)" New conv_gprimd:", conv_gprimd
594 
595    case ("mm2")
596      call wrtout(std_out," Changing the conventional cell: unconventional orthorhombic setting --> conventional","COLL")
597      ! Transformation matrices: unconvential orthorhombic --> conventional orthorhombic. Table 5.1.3.1 ITA page 81.
598      pmat1 = RESHAPE( (/ 0, 1, 0, 1, 0, 0, 0, 0,-1/), (/3,3/) )  ! ( b, a,-c) --> (a,b,c)
599      pmat2 = RESHAPE( (/ 0, 1, 0, 0, 0, 1, 1, 0, 0/), (/3,3/) )  ! ( c, a, b) --> (a,b,c)
600      pmat3 = RESHAPE( (/ 0, 0, 1, 0, 1, 0,-1, 0, 0/), (/3,3/) )  ! (-c, b, a) --> (a,b,c)
601      pmat4 = RESHAPE( (/ 0, 0, 1, 1, 0, 0, 0, 1, 0/), (/3,3/) )  ! ( b, c, a) --> (a,b,c)
602      pmat5 = RESHAPE( (/ 1, 0, 0, 0, 0, 1, 0,-1, 0/), (/3,3/) )  ! ( a,-c, b) --> (a,b,c)
603      conv_gprimd = MATMUL(conv_gprimd,pmat2)
604      !write(std_out,*)" New conv_gprimd:", conv_gprimd
605    case default
606      continue
607    end select
608 
609    ABI_MALLOC(conv_symrec,(3,3,esymm%nsym_gk))
610    conv_symrec = sgk
611 
612    !axes=zero; axes(1,1)=one ; axes(2,2)=one ; axes(3,3)=one
613    !call symrelrot(esymm%nsym_gk,conv_gprimd,axes,conv_symrec,tolsym)
614    call symrelrot(esymm%nsym_gk,Cryst%gprimd,conv_gprimd,conv_symrec,tolsym)
615 
616    ! 3) Reorder indeces such that symmetries are packed in classes.
617    ABI_MALLOC(found,(esymm%nsym_gk))
618    ABI_MALLOC(new_idx,(esymm%nsym_gk))
619    ABI_MALLOC(new_g0,(3,esymm%nsym_gk))
620    new_g0=0; found=0
621 
622    do isym=1,esymm%nsym_gk
623      do jsym=1,esymm%nsym_gk
624        if (ALL(Ptg%sym(:,:,isym) == conv_symrec(:,:,jsym) ))  then
625          found(isym)    = found(isym) + 1
626          new_idx(isym)  = esymm%sgk2symrec(jsym)
627          new_g0(:,isym) = esymm%g0(:,jsym)
628          !EXIT
629        end if
630      end do
631    end do
632    !
633    ! DEBUGGING SECTION
634    !do isym=1,esymm%nsym_gk
635    !  jsym=esymm%sgk2symrec(isym)
636    !  call print_symmetries(1,Cryst%symrec(:,:,jsym),Cryst%tnons(:,jsym),Cryst%symafm(jsym))
637    !  write(std_out,*)esymm%g0(:,isym)
638    !end do
639 
640    if ( Ptg%nsym/=esymm%nsym_gk .or. ANY(found/=1) ) then
641      !write(std_out,*)Ptg%nsym, esymm%nsym_gk
642      !write(std_out,'(a,(i2))')" found = ",found
643      write(std_out,*)" Ptg%sym list, conv_symrec list,  found Ptg% "
644      do isym=1,Ptg%nsym
645        write(std_out,'(a,i2,a,9i2,4x,a,9i2)')" found ",found(isym)," Ptg ",Ptg%sym(:,:,isym),"conv_symrec ",conv_symrec(:,:,isym)
646      end do
647      msg = " sgk and esymm%Ptg are inconsistent. Check tables or source"
648      MSG_WARNING(msg)
649      esymm%err_msg = msg(1:500)
650      esymm%err_status = ESYM_PTG_WRONG_MAPPING
651      esymm%has_chtabs = .FALSE.
652 
653    else ! Reorder symmetries.
654      esymm%sgk2symrec = new_idx
655      esymm%g0 = new_g0
656    end if
657 
658    ABI_FREE(new_idx)
659    ABI_FREE(new_g0)
660    ABI_FREE(found)
661    ABI_FREE(conv_symrec)
662 
663    if (esymm%has_chtabs) then
664      ! Multiply the point group irreps by e^{-ik.\tau} to have the irreps of the little group.
665      ! Store the results in esymm%Ref_irreps so that one can classify the states afterwards.
666      ABI_DT_MALLOC(esymm%Ref_irreps,(esymm%nclass))
667      ABI_MALLOC(phase_mkt,(esymm%nsym_gk))
668 
669      do isym_gk=1,esymm%nsym_gk
670        isym =  esymm%sgk2symrec(isym_gk)
671        mkt = -two_pi * DOT_PRODUCT(esymm%kpt, Cryst%tnons(:,isym))
672        phase_mkt(isym_gk) = CMPLX(DCOS(mkt), DSIN(mkt))
673      end do
674 
675      call copy_irrep(Ptg%Irreps,esymm%Ref_irreps,phase_mkt)
676      ABI_FREE(phase_mkt)
677    end if
678 
679 #if 0
680    ! Herring test requires the evaluation of the expression:
681    !
682    !   sum_{S,\tau} \chi^{k,\alpha} ({S|\tau}^2)
683    !
684    ! where Sk = -k + g0, and \chi is the trace of the \alpha-th
685    ! irreducible representation of the little group of k.
686    ! \chi^{k,\alpha} = e^{-ik.\tau} \chi(\alpha) provided that
687    ! we are not at zone border with a non-symmorphic operation.
688    ! The expression is always real and it can only be equal to \pm Ptg%nsym or zero.
689    ! FIXME this part has to be rewritten from scratch.
690    !if (esymm%err_status/=esymm_NOERROR) then
691    !  write(std_out,*)" Skipping Herring test"
692    !  goto 110
693    !end if
694 
695    if (esymm%can_use_tr) then
696      ABI_MALLOC(her_test,(esymm%nclass))
697 
698      ABI_MALLOC(tr_conv_symrec,(3,3,esymm%nsym_trgk))
699      do isym_tr=1,esymm%nsym_trgk
700        isym = esymm%tr_sgk2symrec(isym_tr)
701        tr_conv_symrec(:,:,isym_tr)=Cryst%symrec(:,:,isym)
702      end do
703 
704      call symrelrot(esymm%nsym_trgk,Cryst%gprimd,conv_gprimd,tr_conv_symrec_tr,tolsym)
705 
706      do isym_tr=1,esymm%nsym_trgk
707        isym = esymm%tr_sgk2symrec(isym_tr)
708        !rot2 = MATMUL(tr_sgk(:,:,isym),tr_sgk(:,:,isym))
709        !tau2 = MATMUL(tr_sgk(:,:,isym),Cryst%tnons(:,isym)) + Cryst%tnons(:,isym)
710 
711        rot2 = MATMUL(tr_conv_symrec(:,:,isym_tr),tr_conv_symrec(:,:,isym_tr))
712        tau2 = MATMUL(tr_conv_symrec(:,:,isym_tr),Cryst%tnons(:,isym)) + Cryst%tnons(:,isym)
713 
714        phase_k = EXP(-j_dpc*two_pi*DOT_PRODUCT(kpoint,tau2))
715        call locate_sym(Ptg,rot2,isym2,icls2)
716 
717        do irp=1,esymm%nclass
718          her_test(irp) = her_test(irp) + phase_k * Ptg%Irreps(irp)%trace(icls2)
719        end do
720      end do
721 
722      ABI_FREE(tr_conv_symrec)
723 
724      ! FIXME
725      ABI_MALLOC(esymm%herring_test,(esymm%nclass))
726 
727      do irp=1,esymm%nclass
728        if ( ABS(her_test(irp) - Ptg%nsym) < tol6 ) then
729          esymm%herring_test(irp) = +1
730        else if ( ABS(her_test(irp)) < tol6 ) then
731          esymm%herring_test(irp) =  0
732        else if ( ABS(her_test(irp) + Ptg%nsym) < tol6 ) then
733          esymm%herring_test(irp) = -1
734        else
735          write(msg,'(a,i2,2a,i0,a,i2)')&
736 &          "Herring test for the irreducible representation number ",irp,ch10,&
737 &          "gave ",esymm%herring_test(irp),", while it should be 0 or +- ",Ptg%nsym
738           MSG_WARNING(msg)
739           esymm%err_msg   =msg
740           esymm%err_status=esymm_HERRING_WRONG_TEST
741        end if
742      end do
743 
744      ABI_FREE(her_test)
745    end if ! can_use_tr
746 #endif
747    !
748    ! Final check
749    !allocate(mtab(esymm%nsym_gk,esymm%nsym_gk))
750    !call mult_table(esymm%nsym_gk,Ptg%sym,mtab)
751 
752    !do isym=1,esymm%nsym_gk
753    !  isym1 = esymm%sgk2symrec(isym)
754    !  do jsym=1,esymm%nsym_gk
755    !    isym2 = esymm%sgk2symrec(jsym)
756    !    rot2 = MATMUL(Cryst%symrec(:,:,isym1),Cryst%symrec(:,:,isym2))
757 
758    !    iprod = mtab(isym,jsym)
759 
760    !    do irp=1,esymm%nclass
761    !       dim_irrep = Ptg%Irreps(irp)%dim
762    !       allocate(mat_test(dim_irrep,dim_irrep))
763    !       mat_test = Ptg%Irreps(irp)%mat(:,:,isym) * Ptg%Irreps(irp)%mat(:,:,jsym)
764    !       !call locate_sym(Ptg,rot2,isym2,icls2)
765    !       write(std_out,*)mat_test - Ptg%Irreps(irp)%mat(:,:,iprod)
766    !       deallocate(mat_test)
767    !    end do
768    !
769    !  end do
770    !end do
771    !
772    !deallocate(mtab)
773  end if
774 
775  ABI_FREE(sgk)
776  ABI_FREE(tr_sgk)
777  ABI_FREE(elements_idx)
778 
779  !% allocate(esymm%irrep2b(0:esymm%nclass))
780  !% call nullify_coeff(esymm%irrep2b)
781  !
782  ! 1) Allocate space for the irreducible representations.
783 
784  ! 2) Try to determine if we are in presence of an accidental degeneracy. Sufficient condition:
785  !    There exists a set of degenerate states whose dimension is greater than the dimension
786  !    of the irreducible representations of the point group. The check can be done only
787  !    if Character tables are available.
788 
789  if (esymm%has_chtabs) then
790    ABI_MALLOC(dim_irreps,(esymm%nclass))
791    dim_irreps = (/(esymm%Ref_irreps(irp)%dim, irp=1,esymm%nclass)/)
792  end if
793 
794  nacc_deg=0
795  ABI_MALLOC(esymm%degs_dim,(esymm%ndegs))
796  ABI_DT_MALLOC(esymm%Calc_irreps,(esymm%ndegs))
797 
798  if (esymm%can_use_tr)  then
799    ABI_DT_MALLOC(esymm%trCalc_irreps,(esymm%ndegs))
800  end if
801 
802  do idg=1,esymm%ndegs
803    dim_degs=esymm%degs_bounds(2,idg)-esymm%degs_bounds(1,idg)+1
804 
805    if (esymm%has_chtabs) then
806      if (ALL(dim_degs /= dim_irreps)) then ! An accidental degeneracy is present.
807        nacc_deg=nacc_deg+1
808      end if
809    end if
810 
811    esymm%degs_dim(idg) = dim_degs
812 
813    call init_irrep(esymm%Calc_irreps(idg),esymm%nsym_gk,dim_degs)
814    if (esymm%can_use_tr) then
815      call init_irrep(esymm%trCalc_irreps(idg),esymm%nsym_trgk,dim_degs)
816    end if
817  end do ! idg
818 
819  if (esymm%has_chtabs) then
820    ABI_FREE(dim_irreps)
821    if (nacc_deg/=0) then
822      write(msg,'(a,i0,a)')" Detected ",nacc_deg," accidental degeneracies."
823      MSG_WARNING(msg)
824      esymm%err_status=ESYM_ACCDEG_ERROR
825      ! TODO this should signal to the caller that we have to decompose the calculated representation.
826      esymm%err_msg =msg(1:500)
827    end if
828  end if
829 
830  DBG_EXIT("COLL")
831 
832 end subroutine esymm_init

m_esymm/esymm_print [ Functions ]

[ Top ] [ m_esymm ] [ Functions ]

NAME

 esymm_print

FUNCTION

INPUTS

OUTPUT

PARENTS

      classify_bands

CHILDREN

      xgeev,xginv,zpotrf,ztrsm

SOURCE

855 subroutine esymm_print(esymm,unit,mode_paral,prtvol)
856 
857 
858 !This section has been created automatically by the script Abilint (TD).
859 !Do not modify the following lines by hand.
860 #undef ABI_FUNC
861 #define ABI_FUNC 'esymm_print'
862 !End of the abilint section
863 
864  implicit none
865 
866 !Arguments ------------------------------------
867 !scalars
868  integer,optional,intent(in) :: prtvol,unit
869  character(len=4),optional,intent(in) :: mode_paral
870  type(esymm_t),intent(in) :: esymm
871 
872 !Local variables-------------------------------
873 !scalars
874  integer :: icl,idg,my_unt,my_prtvol
875  integer :: irr_idx,nstates,nunknown,istart,istop,ii
876  character(len=4) :: my_mode
877  character(len=1000) :: fmt,msg,msg0
878 
879 ! *********************************************************************
880 
881  my_unt   =std_out; if (PRESENT(unit      )) my_unt   =unit
882  my_prtvol=0      ; if (PRESENT(prtvol    )) my_prtvol=prtvol
883  my_mode  ='COLL' ; if (PRESENT(mode_paral)) my_mode  =mode_paral
884 
885  write(fmt,*)'(2a,3f8.4,3a,i4,2a,i3,2a,i2,2a,i2,a,',esymm%nclass,'i2,a)'
886  write(msg,fmt)ch10,&
887 &  ' ===== Character of bands at k-point: ',esymm%kpt,' ===== ',ch10,&
888 &  '   Total number of bands analyzed .................. ',esymm%nbnds,ch10,&
889 &  '   Number of degenerate sets detected .............. ',esymm%ndegs,ch10,&
890 &  '   Number of operations in the little group of k ... ',esymm%nsym_gk,ch10,&
891 &  '   Number of classes (irreps) in the group of k .... ',esymm%nclass,' (',(esymm%nelements(icl),icl=1,esymm%nclass),' )'
892  call wrtout(my_unt,msg,my_mode)
893 
894  if (esymm%nonsymmorphic_at_zoneborder) then
895    call wrtout(my_unt," Non-symmorphic small group at zone border. Character analysis not available ",my_mode)
896  end if
897 
898  if (esymm_failed(esymm)) then
899    write(std_out,'(3a)')"Band classification algorithm failed with the error:",ch10,TRIM(esymm%err_msg)
900    write(msg,'(3a)')"Band classification algorithm failed with the error:",ch10,TRIM(esymm%err_msg)
901    call wrtout(my_unt,msg,my_mode)
902  end if
903 
904  !nunknown=0
905  !do iband=1,esymm%nbnds
906  !  irr_idx = esymm%b2irrep(iband)
907  !  if (irr_idx /= 0) then
908  !    if (     esymm%has_chtabs) irr_name = esymm%Ref_Irreps(irr_idx)%name
909  !    if (.not.esymm%has_chtabs) write(irr_name,'(i0)')irr_idx ! use the index instead of the name.
910  !  else
911  !    irr_name = "???"
912  !    nunknown = nunknown +1
913  !  end if
914  !  write(msg,'(a,i3,2a)')' Band ',iband,' belongs to irrep ',TRIM(irr_name)
915  !  call wrtout(my_unt,msg,my_mode)
916  !end do
917 
918  do irr_idx=1,esymm%nclass
919    nstates = size(esymm%irrep2b(irr_idx)%value)
920    if (esymm%has_chtabs) then
921      write(msg0,'(a,i0,3a)')"  Found ",nstates," states with character ",TRIM(esymm%Ref_irreps(irr_idx)%name),": "
922    else
923      write(msg0,'(2(a,i0),a)')"   Found ",nstates," states with character index ",irr_idx,": "
924    end if
925    do istart=1,nstates,20
926      istop=istart+11; if (istop>nstates) istop=nstates
927      write(msg,'(20(1x,i0))')(esymm%irrep2b(irr_idx)%value(ii), ii=istart,istop)
928      if (istart==1) msg = TRIM(msg0)//TRIM(msg)
929      if (istart/=1) msg = "   "//TRIM(msg)
930      call wrtout(my_unt,msg,my_mode)
931    end do
932  end do
933 
934  nunknown = size(esymm%irrep2b(0)%value)
935  if (nunknown > 0) then
936    write(msg0,'(a,i0,a)')" WARNING: ",nunknown," states have not been classified:"
937    do istart=1,nunknown,20
938      istop=istart+11; if (istop>nunknown) istop=nunknown
939      write(msg,'(20(1x,i0))')(esymm%irrep2b(0)%value(ii), ii=istart,istop)
940      if (istart==1) msg = TRIM(msg0)//TRIM(msg)
941      if (istart/=1) msg = "   "//TRIM(msg)
942      call wrtout(my_unt,msg,my_mode)
943    end do
944  end if
945 
946  if (my_prtvol>0 .or. nunknown>0 .or. .not.esymm%has_chtabs) then ! print the calculated character table.
947    call wrtout(my_unt,ch10//" Calculated character table ",my_mode)
948    !write(fmt,*)'(i2,a,i2,1x,',esymm%nclass,'(a,2f6.3),a)'
949    write(fmt,*)'(i2,a,i2,1x,',esymm%nclass,'(a,2f5.2),a)'
950    do idg=1,esymm%ndegs
951      write(msg,fmt)&
952 &      esymm%degs_bounds(1,idg),'-',esymm%degs_bounds(2,idg),&
953 &      ('|',esymm%Calc_irreps(idg)%trace(esymm%nelements(icl)), icl=1,esymm%nclass),'|'
954      call wrtout(my_unt,msg,my_mode)
955    end do
956  end if
957 
958 end subroutine esymm_print

m_esymm/esymm_symmetrize_mels [ Functions ]

[ Top ] [ m_esymm ] [ Functions ]

NAME

  esymm_symmetrize_mels

FUNCTION

INPUTS

  esymm<esymm_t>

PARENTS

      calc_sigc_me,calc_sigx_me,cohsex_me

CHILDREN

      xgeev,xginv,zpotrf,ztrsm

SOURCE

1433 subroutine esymm_symmetrize_mels(esymm,lbnd,ubnd,in_me,out_me)
1434 
1435 
1436 !This section has been created automatically by the script Abilint (TD).
1437 !Do not modify the following lines by hand.
1438 #undef ABI_FUNC
1439 #define ABI_FUNC 'esymm_symmetrize_mels'
1440 !End of the abilint section
1441 
1442  implicit none
1443 
1444 !Arguments ------------------------------------
1445 !scalars
1446  integer :: lbnd,ubnd
1447  type(esymm_t),target,intent(in) :: esymm
1448 !arrays
1449  complex(dpc),intent(in) :: in_me(2,lbnd:ubnd,lbnd:ubnd)
1450  complex(dpc),intent(out) :: out_me(lbnd:ubnd,lbnd:ubnd)
1451 
1452 !Local variables-------------------------------
1453 !scalars
1454  integer :: idg1,b1_start,b1_stop,irp1
1455  integer :: idg2,b2_start,b2_stop,irp2
1456  integer :: ii,jj,ib,jb,kk,kb,lb,ll
1457  complex(dpc) :: tr_ofd,ofd,dsd,tr_dsd
1458  type(irrep_t),pointer :: Irrep1,Irrep2
1459  type(irrep_t),pointer :: tr_Irrep1,tr_Irrep2
1460 
1461 ! *********************************************************************
1462 
1463  if (esymm_failed(esymm)) then
1464    MSG_ERROR("Symmetrization cannot be performed. You should not be here!")
1465  end if
1466 
1467  do idg1=1,esymm%ndegs  ! First loop over set of degenerate states.
1468    b1_start = esymm%degs_bounds(1,idg1)
1469    b1_stop  = esymm%degs_bounds(2,idg1)
1470 
1471    !if (b1_stop<lbnd .or. b2_start >ubnd) then
1472    !  MSG_ERROR("Wrong band indeces, check esymm initialization")
1473    !end if
1474 
1475    Irrep1 => esymm%Calc_irreps(idg1)
1476    if (esymm%can_use_tr) tr_Irrep1 => esymm%trCalc_irreps(idg1)
1477    irp1 = esymm%b2irrep(b1_start)
1478 
1479    do idg2=1,esymm%ndegs ! Second loop over set of degenerate states.
1480      !write(std_out,*)" ==> Symmetrizing degenerate set ",idg1,idg2
1481      b2_start = esymm%degs_bounds(1,idg2)
1482      b2_stop  = esymm%degs_bounds(2,idg2)
1483      irp2 = esymm%b2irrep(b2_start)
1484 
1485      if (irp1/=irp2 .or. idg1==idg2) CYCLE  ! Skip diago elements or elements belonging to different irreps.
1486 
1487      Irrep2 => esymm%Calc_irreps(idg2)
1488      if (esymm%can_use_tr) tr_Irrep2 => esymm%trCalc_irreps(idg2)
1489      !
1490      ! Symmetrize the off-diagonal matrix elements.
1491      ! summing over kk and ll. ii and jj are the indeces of the bands that are symmetrized
1492      do ii=1,b1_stop-b1_start+1
1493        ib= ii+b1_start-1
1494        do jj=1,b2_stop-b2_start+1
1495          jb= jj+b2_start-1
1496          !write(std_out,*)" ====> Symmetrizing ",ib,jb
1497 
1498          ofd= czero; tr_ofd=czero
1499          do kk=1,b1_stop-b1_start+1
1500            kb= kk+b1_start-1
1501            do ll=1,b2_stop-b2_start+1
1502              lb= ll+b2_start-1
1503              dsd = sum_irreps(Irrep1,Irrep2,kk,ii,ll,jj)
1504              ofd = ofd + dsd * in_me(1,kb,lb)
1505              if (esymm%can_use_tr) then
1506                tr_dsd = sum_irreps(tr_Irrep1,tr_Irrep2,kk,jj,ll,ii) ! Exchange of band indeces.
1507                tr_ofd = tr_ofd + tr_dsd * in_me(2,kb,lb)            ! Contribution obtained from TR.
1508              end if
1509            end do
1510          end do
1511 
1512          out_me(ib,jb)= ofd/esymm%nsym_gk
1513          if (esymm%can_use_tr .and. esymm%nsym_trgk>0) out_me(ib,jb)= out_me(ib,jb) + tr_ofd/esymm%nsym_trgk
1514        end do
1515      end do
1516    end do
1517  end do
1518 
1519 end subroutine esymm_symmetrize_mels

m_esymm/esymm_t [ Types ]

[ Top ] [ m_esymm ] [ Types ]

NAME

  esymm_t

FUNCTION

  Dataype gathering data and tables needed to analize the symmetries
  of electronic states at a given k-point via Group Theory.

SOURCE

 73  type,public :: esymm_t
 74 
 75   integer :: nspinor
 76   ! Number of spinorial components.
 77 
 78   integer :: first_ib
 79   ! Index of the first treated band.
 80 
 81   integer :: nbnds
 82   ! Number of bands for this k-point and spin.
 83 
 84   integer :: nclass
 85   ! The number of classes in the group of k.
 86 
 87   integer :: nsym_gk
 88   ! Number of symmetries in the group of k. Namely that the set of symmetries such that Sk = k +G0.
 89 
 90   integer :: nsym_trgk
 91   ! Number of symmetries in the extended group of k. Namely that the set of symmetries such that -Sk = k + G0.
 92 
 93   integer :: err_status = ESYM_NOERROR
 94   ! Flag signaling if the classification algorithm succeed or not.
 95 
 96   real(dp) :: tol_deg
 97   ! Energy tolerance below which two states are considered degenerate.
 98 
 99   logical :: can_use_tr
100   ! .TRUE. if time-reversal can be used
101 
102   logical :: only_trace
103   ! if .TRUE. only the trace of a single matrix per class is calculated
104   ! this is the standard way used to analyze bands symmetries. If .FALSE.
105   ! the full matrices of the irreducible representations are calculated and stored
106 
107   logical :: has_spatial_inv
108   ! .TRUE. if the inversion belongs to the space group
109 
110   logical :: nonsymmorphic_at_zoneborder
111   ! if .TRUE. analysis cannot be performed since kpt is
112   ! at border zone and non-zero fractional translations are present in the space group
113 
114   logical :: has_chtabs
115   ! True if Ref_irreps and character tables are available (tables are initialized either
116   ! from point group irreps or from an external database downloaded from the Bilbao server)
117 
118   real(dp) :: kpt(3)
119   ! The crystalline momentum of the wavefunctions in reduced coordinates.
120 
121   character(len=500) :: err_msg="None"
122 
123   integer,allocatable :: g0(:,:)
124   ! g0(3,nsym_gk)
125   ! The umklapp g0 vector associated to each little group operation.
126 
127   integer,allocatable :: tr_g0(:,:)
128   ! tr_g0(3,nsym_trgk)
129   ! The umklapp g0 vector associated to each little group operation.
130 
131   integer :: ndegs
132   ! Number of degenerate states.
133 
134   integer,allocatable :: nelements(:)
135   ! nelements(nclass)
136   ! Number of symmetry operations in each class.
137 
138   integer,allocatable :: sgk2symrec(:)
139   ! sgk2symrec(nsym_gk)
140   ! Mapping between the symmetries of the group of k and the symrec(l) array.
141   ! The symmetries of the little group are always packed in classes to facilitate
142   ! the calculation of the character of the irrep. Abinit symmetries are randomly ordered.
143 
144   integer,allocatable :: tr_sgk2symrec(:)
145   ! trsgk2symrec(nsym_trgk)
146   ! Mapping between the symmetries of the group of k and the symrec(l) array.
147   ! The symmetries of the little group are always packed in classes to facilitate
148   ! the calculation of the character of the irrep. Abinit symmetries are randomly ordered.
149 
150   integer,allocatable :: herring_test(:)
151   ! herring_test(nclass)
152   ! The result of Herring test for each irreducible representantion of the group of k.
153   ! Possible values are: +1, 0, -1
154 
155   integer,allocatable :: b2irrep(:)
156   ! b2irrep(nbnds)
157   ! For each band, it gives the index of the irreducible representation in Ref_irreps.
158 
159   type(coeffi1_type),allocatable :: irrep2b(:)
160   ! irrep2b(0:nclass)%value(:)
161   ! Ragged arrays with the mapping between the set of irreducible representation and the band indices.
162   ! irrep2b(irp)%value(:) gives the indeces of the states belonging to irrep irp, irp=1,nclass
163   ! irrep2b(0)%value(:) stores the indeces of the states that have not been classified due to
164   !   the presence of an accidental degeneracy.
165 
166   integer,allocatable :: degs_bounds(:,:)
167   ! degs_bounds(2,ndegs)
168   !   degs_bounds(1,idg)= first band index of the degenerate set idg=1,ndegs
169   !   degs_bounds(2,idg)= final band index of the degenerate set idg=1,ndegs
170 
171   integer,allocatable :: degs_dim(:)
172   ! degs_dim(ndegs)
173   ! Number of states in each degenerate subspace. Cannot be larger that nclass provided
174   ! that no accidental degeneracy occurs.
175 
176   !% integer,allocatable :: class_ids(:,:)
177   ! class_ids(2,nclass)
178   ! (1,icl) = index of the first symmetry of class icl
179   ! (2,icl) = index of the last symmetry of class icl
180   ! Note that symmetries in sym are packed in classes.
181 
182   type(irrep_t),allocatable :: Calc_irreps(:)
183   ! Calc_irreps(ndegs)
184   !  The representations of the little group of k calculated from the wavefunctions. <\phi_nk|R_t|\phi_mk>
185   !  where R_t belong to the little group of k.
186   !  They represent an unitary irreducible representation provided that no accidental degeneracy occurs.
187 
188   type(irrep_t),allocatable :: trCalc_irreps(:)
189   ! trCalc_irreps(ndegs)
190   !  The representations of the little group of k calculated from the wavefunctions. <\phi_nk|R_t|\phi_mk>
191   !  where R_t belong to the little group of k.
192   !  They represent an unitary irreducible representation provided that no accidental degeneracy occurs.
193 
194   type(irrep_t),allocatable :: Ref_irreps(:)
195   ! Irreps(nclass)
196   !   Reference irreducible representations of the group of k derived from the point group
197   !   or from the external database downloaded from the Bilbao web site.
198 
199  end type esymm_t
200 
201  public :: esymm_init             ! Initialize the object
202  public :: esymm_print            ! Print info
203  public :: esymm_free             ! Free memory
204  public :: esymm_finalize         ! Finalize the object
205  public :: esymm_symmetrize_mels  ! Symmetrize given matrix elements
206  public :: esymm_failed           ! True if symmetry analysis failed.

m_esymm/polish_irreps [ Functions ]

[ Top ] [ m_esymm ] [ Functions ]

NAME

  polish_irreps

FUNCTION

INPUTS

PARENTS

      m_esymm

CHILDREN

      xgeev,xginv,zpotrf,ztrsm

SOURCE

1578 subroutine polish_irreps(Irreps)
1579 
1580 
1581 !This section has been created automatically by the script Abilint (TD).
1582 !Do not modify the following lines by hand.
1583 #undef ABI_FUNC
1584 #define ABI_FUNC 'polish_irreps'
1585 !End of the abilint section
1586 
1587  implicit none
1588 
1589 !Arguments ------------------------------------
1590 !scalars
1591  type(irrep_t),intent(inout) :: Irreps(:)
1592 
1593 !Local variables-------------------------------
1594 !scalars
1595  integer,parameter :: ldvl1=1,ldvr1=1
1596  integer :: irp,sym,dim,ldvr,ii,ivec,jvec,info
1597  !character(len=500) :: msg
1598 !arrays
1599  complex(dpc),allocatable :: vl(:,:),vr(:,:),vrm1(:,:),overlap(:,:)
1600  complex(dpc),allocatable :: cmat(:,:),eigval(:)
1601 
1602 ! *********************************************************************
1603 
1604  ! Eigen decomposition: A = V D V^{-1}.
1605  do irp=1,SIZE(Irreps)
1606    dim = Irreps(irp)%dim
1607    ABI_MALLOC(cmat,(dim,dim))
1608    ABI_MALLOC(eigval,(dim))
1609    ldvr=dim
1610    ABI_MALLOC(vl,(ldvl1,dim))
1611    ABI_MALLOC(vr,(ldvr,dim))
1612    ABI_MALLOC(vrm1,(dim,dim))
1613    ABI_MALLOC(overlap,(dim,dim))
1614    do sym=1,Irreps(irp)%nsym
1615      cmat = Irreps(irp)%mat(:,:,sym)
1616      call xgeev("No vectors","Vectors",dim,cmat,dim,eigval,vl,ldvl1,vr,ldvr)
1617      !
1618      ! Orthogonalize the eigenvectors using Cholesky orthogonalization.
1619      do jvec=1,dim
1620        do ivec=1,jvec
1621          overlap(ivec,jvec) = DOT_PRODUCT(vr(:,ivec),vr(:,jvec))
1622        end do
1623      end do
1624      !
1625      ! 2) Cholesky factorization: overlap = U^H U with U upper triangle matrix.
1626      call ZPOTRF('U',dim,overlap,dim,info)
1627      ABI_CHECK(info == 0, sjoin('ZPOTRF returned info=', itoa(info)))
1628 
1629      ! 3) Solve X U = Vr, on exit the Vr treated by this node is orthonormalized.
1630      call ZTRSM('R','U','N','N',dim,dim,cone,overlap,dim,vr,dim)
1631 
1632      !write(std_out,*)"After ortho",MATMUL(TRANSPOSE(CONJG(vr)),vr)
1633 
1634      vrm1 = vr
1635      call xginv(vrm1,dim)
1636      do ii=1,dim
1637        eigval(ii) = eigval(ii)/ABS(eigval(ii)) ! Rescale the eigevalues.
1638        vrm1(ii,:) =  eigval(ii) * vrm1(ii,:)
1639      end do
1640      Irreps(irp)%mat(:,:,sym) = MATMUL(vr,vrm1)
1641      Irreps(irp)%trace(sym) = get_trace(Irreps(irp)%mat(:,:,sym))
1642    end do
1643    ABI_FREE(cmat)
1644    ABI_FREE(eigval)
1645    ABI_FREE(vl)
1646    ABI_FREE(vr)
1647    ABI_FREE(vrm1)
1648    ABI_FREE(overlap)
1649  end do
1650 
1651 end subroutine polish_irreps

m_esymm/which_irrep [ Functions ]

[ Top ] [ m_esymm ] [ Functions ]

NAME

  m_esymm

FUNCTION

  Return the index of the irreducible representation with character charact. 0 if not found.

INPUTS

  esymm<esymm_t>
  trace(%nsym_gk)=The trace of the representation to be compared with the internal database (if present).
  tolerr=Absolute error on the character.

PARENTS

SOURCE

1376 function which_irrep(esymm,trace,tolerr)
1377 
1378 
1379 !This section has been created automatically by the script Abilint (TD).
1380 !Do not modify the following lines by hand.
1381 #undef ABI_FUNC
1382 #define ABI_FUNC 'which_irrep'
1383 !End of the abilint section
1384 
1385  implicit none
1386 
1387 !Arguments ------------------------------------
1388 !scalars
1389  integer :: which_irrep
1390  real(dp),intent(in) :: tolerr
1391  type(esymm_t),intent(in) :: esymm
1392 !arrays
1393  complex(dpc),intent(in) :: trace(esymm%nsym_gk)
1394 
1395 !Local variables-------------------------------
1396 !scalars
1397  integer :: irp
1398 
1399 ! *********************************************************************
1400 
1401  which_irrep = 0
1402  if (esymm%has_chtabs) then ! Symmetry analysis can be performed.
1403    do irp=1,esymm%nclass
1404      if ( ALL( ABS(esymm%Ref_irreps(irp)%trace(:) - trace(:)) < tolerr)) then
1405        which_irrep = irp; EXIT
1406      end if
1407    end do
1408  end if
1409 
1410 end function which_irrep