TABLE OF CONTENTS
- ABINIT/m_esymm
- m_esymm/esymm_failed
- m_esymm/esymm_finalize
- m_esymm/esymm_free_0D
- m_esymm/esymm_free_2D
- m_esymm/esymm_init
- m_esymm/esymm_print
- m_esymm/esymm_symmetrize_mels
- m_esymm/esymm_t
- m_esymm/polish_irreps
- m_esymm/which_irrep
ABINIT/m_esymm [ 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-2024 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 .
SOURCE
17 #if defined HAVE_CONFIG_H 18 #include "config.h" 19 #endif 20 21 #include "abi_common.h" 22 23 MODULE m_esymm 24 25 use defs_basis 26 use m_abicore 27 use m_errors 28 29 use m_io_tools, only : file_exists 30 use m_symtk, only : matr3inv, sg_multable, symrelrot, littlegroup_q 31 use m_symfind, only : symbrav 32 use m_fstrings, only : int2char10, itoa, sjoin 33 use m_numeric_tools, only : print_arr, set2unit, get_trace 34 use m_hide_lapack, only : xgeev, xginv 35 use m_crystal, only : crystal_t 36 use m_defs_ptgroups, only : point_group_t, irrep_t 37 use m_ptgroups, only : get_classes, point_group_init, irrep_free,& 38 & copy_irrep, init_irrep, mult_table, sum_irreps 39 40 implicit none 41 42 private
m_esymm/esymm_failed [ Functions ]
[ Top ] [ m_esymm ] [ Functions ]
NAME
esymm_failed
FUNCTION
INPUTS
esymm<esymm_t>
SOURCE
1404 function esymm_failed(esymm) 1405 1406 !Arguments ------------------------------------ 1407 !scalars 1408 logical :: esymm_failed 1409 type(esymm_t),intent(in) :: esymm 1410 1411 ! ********************************************************************* 1412 1413 esymm_failed = (esymm%err_status /= ESYM_NOERROR) 1414 1415 end function esymm_failed
m_esymm/esymm_finalize [ Functions ]
[ Top ] [ m_esymm ] [ Functions ]
NAME
esymm_finalize
FUNCTION
INPUTS
OUTPUT
SOURCE
1017 subroutine esymm_finalize(esymm,prtvol) 1018 1019 !Arguments ------------------------------------ 1020 !scalars 1021 integer,intent(in) :: prtvol 1022 type(esymm_t),target,intent(inout) :: esymm 1023 1024 !Local variables------------------------------- 1025 integer :: idg,ib1,ib2,idx,nunknown,dg_dim 1026 integer :: try,irep,nitems,nseen,isn 1027 integer :: isym,idg1,idg2,dim_mat,irr_idx2,irr_idx1 1028 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. 1029 !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. 1030 !real(dp),parameter :: TOL_TRACE=tol3,TOL_ORTHO=tol3,TOL_UNITARY=tol3 ! Large tolerance is needed to avoid problems. 1031 real(dp) :: uerr,max_err 1032 complex(dpc) :: ctest 1033 logical :: isnew 1034 character(len=500) :: msg 1035 !arrays 1036 integer,allocatable :: dims_seen(:) 1037 complex(dpc),allocatable :: traces_seen(:,:) 1038 complex(dpc),pointer :: trace(:) 1039 complex(dpc),pointer :: calc_mat(:,:),trace1(:),trace2(:) 1040 complex(dpc),allocatable :: cidentity(:,:) 1041 1042 ! ************************************************************************* 1043 1044 !@esymm_t 1045 1046 ! Each band is initialized as "Unknown". 1047 esymm%b2irrep = 0 1048 1049 ! Force the matrices to be unitary. 1050 call polish_irreps(esymm%Calc_irreps) 1051 1052 if (.not.esymm%has_chtabs) then 1053 1054 write(msg,'(5a)')& 1055 & "Reference character table not available. ",ch10,& 1056 & "Symmetry analysis not available. Using heuristic method to classify the states.",ch10,& 1057 & "It might not work, especially if accidental degeneracies are present." 1058 ABI_WARNING(msg) 1059 ! 1060 ! The simplest thing we can do here is using the calculated matrices to get the 1061 ! character and comparing the results hoping everything is OK. 1062 ABI_MALLOC(traces_seen,(esymm%nsym_gk,esymm%ndegs)) 1063 ABI_MALLOC(dims_seen,(esymm%ndegs)) 1064 1065 traces_seen=czero; nseen=1 1066 traces_seen(:,1) = esymm%Calc_irreps(1)%trace 1067 dims_seen(1) = esymm%Calc_irreps(1)%dim 1068 1069 do idg=2,esymm%ndegs 1070 dg_dim = esymm%Calc_irreps(idg)%dim 1071 trace => esymm%Calc_irreps(idg)%trace 1072 isnew=.TRUE. 1073 do isn=1,nseen 1074 if (ALL (ABS(trace - traces_seen(:,isn)) < TOL_TRACE) ) then 1075 isnew=.FALSE.; EXIT 1076 end if 1077 end do 1078 1079 if (isnew) then 1080 nseen = nseen+1 1081 traces_seen(:,nseen) = trace 1082 dims_seen(nseen) = dg_dim 1083 end if 1084 end do 1085 1086 if (nseen>esymm%nclass) then 1087 write(msg,'(3a)')& 1088 & "The number of different calculated traces is found to be greater than nclasses!",ch10,& 1089 & "Heuristic method clearly failed. Symmetry analysis cannot be performed." 1090 ABI_WARNING(msg) 1091 esymm%err_status = ESYM_HEUR_WRONG_NCLASSES 1092 esymm%err_msg = msg 1093 1094 do isn=1,nseen 1095 write(msg,'(a,i0)')" Representation: ",isn 1096 call wrtout(std_out,msg,"COLL") 1097 call print_arr(traces_seen(:,isn),max_r=esymm%nsym_gk,unit=std_out,mode_paral="COLL") 1098 end do 1099 1100 else ! It seems that the Heuristic method succeeded. 1101 do idg=1,esymm%ndegs 1102 ib1=esymm%degs_bounds(1,idg) 1103 ib2=esymm%degs_bounds(2,idg) 1104 trace => esymm%Calc_irreps(idg)%trace 1105 do isn=1,nseen 1106 if (ALL (ABS(trace - traces_seen(:,isn)) < TOL_TRACE) ) then 1107 esymm%b2irrep(ib1:ib2)=isn 1108 if (esymm%Calc_irreps(idg)%dim /= dims_seen(isn)) then 1109 write(msg,'(3a)')& 1110 & "Found two set of degenerate states with same character but different dimension!",ch10,& 1111 & "heuristic method clearly failed. Symmetry analysis cannot be performed." 1112 ABI_ERROR(msg) 1113 esymm%err_status = ESYM_HEUR_WRONG_DIMS 1114 esymm%err_msg = msg 1115 end if 1116 EXIT 1117 end if 1118 end do 1119 end do 1120 end if 1121 1122 ABI_FREE(traces_seen) 1123 ABI_FREE(dims_seen) 1124 1125 else 1126 ! 1127 ! * Search in the lookup table definining the irreducible representation 1128 nunknown = 0 1129 do idg=1,esymm%ndegs 1130 1131 ib1=esymm%degs_bounds(1,idg) 1132 ib2=esymm%degs_bounds(2,idg) 1133 trace => esymm%Calc_irreps(idg)%trace 1134 1135 try = which_irrep(esymm, trace, tol3) 1136 if (try==0) try = which_irrep(esymm, trace, 0.1_dp) ! try again with increased tolerance. 1137 if (try/=0) then 1138 esymm%b2irrep(ib1:ib2)=try 1139 else 1140 esymm%b2irrep(ib1:ib2)=0 1141 nunknown = nunknown + (ib2-ib1+1) 1142 end if 1143 end do 1144 end if 1145 ! 1146 ! %irrep2b(0)) gives the indices of the states that have not been classified. 1147 ABI_MALLOC(esymm%irrep2b,(0:esymm%nclass)) 1148 1149 !write(std_out,*)"b2irrep",esymm%b2irrep 1150 1151 do irep=0,esymm%nclass 1152 nitems = COUNT(esymm%b2irrep==irep) 1153 ABI_MALLOC(esymm%irrep2b(irep)%value,(nitems)) 1154 idx=0 1155 do ib1=1,esymm%nbnds 1156 if (esymm%b2irrep(ib1) == irep) then 1157 idx = idx + 1 1158 esymm%irrep2b(irep)%value(idx) = ib1 1159 end if 1160 end do 1161 end do 1162 1163 if (size(esymm%irrep2b(0)%value) /= 0) then 1164 write(msg,'(a,i0,a)')" Band classification algorithm was not able to classify ",size(esymm%irrep2b(0)%value)," states." 1165 ABI_WARNING(msg) 1166 esymm%err_status = ESYM_CLASSIFICATION_ERROR 1167 esymm%err_msg = msg 1168 end if 1169 ! 1170 ! ============================================================== 1171 ! ==== Test basic properties of irreducible representations ==== 1172 ! ============================================================== 1173 1174 if (.not.esymm_failed(esymm)) then 1175 ! 1176 ! 1) \sum_R \chi^*_a(R)\chi_b(R)= N_R \delta_{ab} 1177 ! 1178 !call wrtout(std_out," \sum_R \chi^*_a(R)\chi_b(R) = N_R \delta_{ab} ","COLL") 1179 max_err=zero 1180 do idg2=1,esymm%ndegs 1181 trace2 => esymm%Calc_irreps(idg2)%trace(1:esymm%nsym_gk) 1182 ib2 = esymm%degs_bounds(1,idg2) 1183 irr_idx2 = esymm%b2irrep(ib2) 1184 if (irr_idx2 == 0) CYCLE 1185 1186 do idg1=1,idg2 1187 trace1 => esymm%Calc_irreps(idg1)%trace(1:esymm%nsym_gk) 1188 ib1 = esymm%degs_bounds(1,idg1) 1189 irr_idx1 = esymm%b2irrep(ib1) 1190 if (irr_idx1 == 0) CYCLE 1191 ctest=DOT_PRODUCT(trace1,trace2)/esymm%nsym_gk 1192 if (irr_idx1==irr_idx2) ctest=ctest-one 1193 max_err = MAX(max_err,ABS(ctest)) 1194 if (.FALSE..and.ABS(ctest)>tol3) then 1195 write(msg,'(a,4i3,2es16.8)')& 1196 & ' WARNING: should be delta_ij: cx1 cx2, irr1, irr2, ctest: ',idg1,idg2,irr_idx1,irr_idx2,ctest 1197 call wrtout(std_out,msg,"COLL") 1198 end if 1199 end do 1200 end do 1201 1202 if (max_err>TOL_ORTHO) then 1203 write(msg,'(a,es10.2)')" Too large maximum error on \sum_R \chi^*_a(R)\chi_b(R) = N_R \delta_{ab}: ",max_err 1204 ABI_WARNING(msg) 1205 esymm%err_status = ESYM_ORTHO_ERROR 1206 esymm%err_msg = msg 1207 else 1208 write(msg,'(a,es10.2)')" maximum error on \sum_R \chi^*_a(R)\chi_b(R) = N_R \delta_{ab}: ",max_err 1209 call wrtout(std_out,msg,"COLL") 1210 end if 1211 1212 if (.not.esymm%only_trace) then 1213 !call wrtout(std_out," **** Testing the unitary of the calculated irreps ****",my_mode) 1214 max_err=zero 1215 do idg1=1,esymm%ndegs 1216 ib1 = esymm%degs_bounds(1,idg1) 1217 irr_idx1 = esymm%b2irrep(ib1) 1218 if (irr_idx1 == 0) CYCLE 1219 1220 do isym=1,esymm%nsym_gk 1221 calc_mat => esymm%Calc_irreps(idg1)%mat(:,:,isym) 1222 dim_mat = esymm%Calc_irreps(idg1)%dim 1223 ABI_MALLOC(cidentity,(dim_mat,dim_mat)) 1224 call set2unit(cidentity) 1225 uerr = MAXVAL( ABS(MATMUL(calc_mat,TRANSPOSE(DCONJG(calc_mat))) - cidentity) ) 1226 max_err = MAX(max_err,uerr) 1227 ABI_FREE(cidentity) 1228 if (.FALSE..and.prtvol>=10) then 1229 write(std_out,'(a,i3,a,i2,a,es16.8,a)')& 1230 & " === idg: ",idg1,", isym: ",isym,", Error on U^* U = 1: ",uerr," ===" 1231 call print_arr(calc_mat,dim_mat,dim_mat,unit=std_out,mode_paral="COLL") 1232 end if 1233 end do 1234 end do 1235 1236 if (max_err>TOL_UNITARY) then 1237 write(msg,'(a,es10.2)')" Too large maximum error on the unitary of representions matrices: ",max_err 1238 ABI_WARNING(msg) 1239 esymm%err_msg = msg 1240 esymm%err_status = ESYM_UNITARY_ERROR 1241 else 1242 write(msg,'(a,es10.2)')" maximum error on the unitary of representions matrices: ",max_err 1243 call wrtout(std_out,msg,"COLL") 1244 end if 1245 1246 end if 1247 1248 end if 1249 1250 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)
SOURCE
937 subroutine esymm_free_0D(esymm) 938 939 !Arguments ------------------------------------ 940 !scalars 941 type(esymm_t),intent(inout) :: esymm 942 943 !Local variables ------------------------------ 944 !scalars 945 integer :: ii 946 ! ************************************************************************* 947 948 !@esymm_t 949 ABI_SFREE(esymm%g0) 950 ABI_SFREE(esymm%tr_g0) 951 ABI_SFREE(esymm%nelements) 952 ABI_SFREE(esymm%sgk2symrec) 953 ABI_SFREE(esymm%tr_sgk2symrec) 954 ABI_SFREE(esymm%herring_test) 955 ABI_SFREE(esymm%b2irrep) 956 ABI_SFREE(esymm%degs_bounds) 957 ABI_SFREE(esymm%degs_dim) 958 959 if (allocated(esymm%irrep2b)) then 960 do ii=LBOUND(esymm%irrep2b,DIM=1),UBOUND(esymm%irrep2b,DIM=1) 961 ABI_FREE(esymm%irrep2b(ii)%value) 962 end do 963 ABI_FREE(esymm%irrep2b) 964 end if 965 966 if (allocated(esymm%Calc_irreps)) call irrep_free(esymm%Calc_irreps) 967 if (allocated(esymm%trCalc_irreps)) call irrep_free(esymm%trCalc_irreps) 968 if (allocated(esymm%Ref_irreps)) call irrep_free(esymm%Ref_irreps) 969 970 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)
SOURCE
984 subroutine esymm_free_2D(esymm) 985 986 !Arguments ------------------------------------ 987 !scalars 988 type(esymm_t),intent(inout) :: esymm(:,:) 989 990 !Local variables ------------------------------ 991 integer :: id1,id2 992 ! ************************************************************************* 993 994 do id2=1,SIZE(esymm,DIM=2) 995 do id1=1,SIZE(esymm,DIM=1) 996 call esymm_free_0D(esymm(id1,id2)) 997 end do 998 end do 999 1000 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.
SOURCE
246 subroutine esymm_init(esymm,kpt_in,Cryst,only_trace,nspinor,first_ib,nbnds,EDIFF_TOL,ene_k,tolsym) 247 248 !Arguments ------------------------------------ 249 !scalars 250 integer,intent(in) :: nbnds,nspinor,first_ib 251 real(dp),intent(in) :: EDIFF_TOL,tolsym 252 logical,intent(in) :: only_trace 253 type(crystal_t),intent(in) :: Cryst 254 type(esymm_t),intent(out) :: esymm 255 !arrays 256 real(dp),intent(in) :: ene_k(nbnds),kpt_in(3) 257 258 !Local variables------------------------------- 259 !scalars 260 integer :: dim_degs,iband,idg,irp,nacc_deg,isym_gk,grp_ierr 261 integer :: nsym_fm,idx_fm,idx_gk,idx_trgk,isym,jsym,dummy_timrev !,iholohedry 262 integer :: iel,icls,msym,iord !isym1,!iprod,dim_irrep,icls2, isym2,isym_tr, 263 integer :: spgroup,chkprim !,ptgroupma 264 real(dp) :: mkt 265 !complex(dpc) :: phase_k 266 character(len=5) :: ptgroup,ptgroup_name 267 character(len=10) :: spgroup_str 268 character(len=1000) :: msg 269 character(len=fnlen) :: lgroup_fname 270 !arrays 271 integer :: inversion(3,3) 272 integer,allocatable :: degs_bounds(:,:),dim_irreps(:) 273 integer :: bravais(11),sym_axis(3) 274 real(dp) :: pmat1(3,3),pmat2(3,3),pmat3(3,3),pmat4(3,3),pmat5(3,3),pmat6(3,3) 275 !real(dp) :: genafm(3) 276 !integer :: rot2(3,3) 277 !integer,allocatable :: mtab(:,:) 278 integer,allocatable :: elements_idx(:,:),tmp_nelements(:) 279 integer,allocatable :: found(:),symrec_fm(:,:,:),fm2symrec(:) 280 integer,allocatable :: ksym_table(:,:,:),sgk(:,:,:),tr_sgk(:,:,:),dum_symafm(:) 281 integer,allocatable :: new_idx(:),new_g0(:,:),tmp_symrec(:,:,:),conv_symrec(:,:,:) !,tr_conv_symrec(:,:,:) 282 integer,allocatable :: dummy_symafm(:) 283 real(dp) :: conv_gprimd(3,3),axes(3,3) !,tau2(3) 284 !complex(dpc),allocatable :: her_test(:) !,mat_test(:,:) 285 complex(dpc),allocatable :: phase_mkt(:) 286 type(point_group_t) :: Ptg 287 288 ! ************************************************************************* 289 290 DBG_ENTER("COLL") 291 292 !@esymm_t 293 esymm%err_status= ESYM_NOERROR 294 inversion=RESHAPE((/-1,0,0,0,-1,0,0,0,-1/),(/3,3/)) 295 ! 296 ! ==================================== 297 ! ==== Initialize basic variables ==== 298 ! ==================================== 299 esymm%nspinor = nspinor 300 esymm%first_ib = first_ib 301 esymm%nbnds = nbnds 302 esymm%only_trace = only_trace 303 esymm%tol_deg = EDIFF_TOL 304 esymm%has_spatial_inv= (cryst%idx_spatial_inversion() /= 0) 305 esymm%can_use_tr = .TRUE. !TODO this should be input 306 esymm%has_chtabs = .FALSE. 307 esymm%kpt = kpt_in(:) 308 esymm%nonsymmorphic_at_zoneborder=.FALSE. 309 ! 310 ! =============================== 311 ! === Locate degenerate_bands === 312 ! =============================== 313 esymm%ndegs=1 314 315 ABI_MALLOC(degs_bounds,(2,nbnds)) 316 degs_bounds=0; degs_bounds(1,1)=1 317 318 do iband=2,nbnds 319 if (ABS(ene_k(iband)-ene_k(iband-1))>EDIFF_TOL) then 320 degs_bounds(2,esymm%ndegs) = iband-1 + (first_ib-1) 321 esymm%ndegs=esymm%ndegs+1 322 degs_bounds(1,esymm%ndegs) = iband + (first_ib-1) 323 end if 324 end do 325 degs_bounds(2,esymm%ndegs)=nbnds + (first_ib-1) 326 327 ABI_MALLOC(esymm%degs_bounds,(2,esymm%ndegs)) 328 esymm%degs_bounds = degs_bounds(:,1:esymm%ndegs) 329 ABI_FREE(degs_bounds) 330 ! 331 ! Each band is initialized as "Unknown". 332 ABI_MALLOC(esymm%b2irrep,(esymm%nbnds)) 333 esymm%b2irrep = 0 334 ! 335 ! ================================== 336 ! ==== Find the group of kpt_in ==== 337 ! ================================== 338 ! The small point group is the subset of symrec such that $ S q = q + g0 $ 339 ! Symmetries are packed in classes. 340 ! For the time being, AFM symmetries are not treated. 341 342 write(msg,'(a,3(1x,f7.4))')" Finding the little group of k-point: ",esymm%kpt 343 call wrtout(std_out,msg,"COLL") 344 345 ! Only FM symmetries are used. 346 nsym_fm = COUNT(Cryst%symafm==1) 347 348 if (nsym_fm /= Cryst%nsym) then 349 write(msg,'(4a)')ch10,& 350 & "Band classification in terms of magnetic space groups not coded! ",ch10,& 351 & "Only the ferromagnetic subgroup will be used " 352 ABI_COMMENT(msg) 353 end if 354 355 ABI_MALLOC(symrec_fm,(3,3,nsym_fm)) 356 ABI_MALLOC(fm2symrec,(nsym_fm)) 357 358 idx_fm = 0 359 do isym=1,Cryst%nsym 360 if (Cryst%symafm(isym) == 1) then 361 idx_fm = idx_fm+1 362 symrec_fm(:,:,idx_fm) = Cryst%symrec(:,:,isym) 363 fm2symrec(idx_fm) = isym 364 end if 365 end do 366 367 ! Find symmetries that preserve k. 368 ABI_MALLOC(ksym_table,(4,2,nsym_fm)) 369 ABI_MALLOC(dummy_symafm,(nsym_fm)) 370 371 dummy_symafm = 1 372 call littlegroup_q(nsym_fm,esymm%kpt,ksym_table,symrec_fm,dummy_symafm,dummy_timrev,prtvol=0) 373 374 esymm%nsym_gk =COUNT(ksym_table(4,1,:)==1) ! # S such that S k = k +G0 375 376 esymm%nsym_trgk=0 377 if (esymm%can_use_tr) esymm%nsym_trgk=COUNT(ksym_table(4,2,:)==1) ! # S such that -S k = k +G0 378 379 ! Allocate workspace arrays. 380 ABI_MALLOC(sgk,(3,3,esymm%nsym_gk)) 381 ABI_MALLOC(tr_sgk,(3,3,esymm%nsym_trgk)) 382 383 ! Allocate mapping little-group --> symrec and table for umklapps. 384 ABI_MALLOC(esymm%sgk2symrec,(esymm%nsym_gk)) 385 ABI_MALLOC(esymm%g0,(3,esymm%nsym_gk)) 386 ABI_MALLOC(esymm%tr_sgk2symrec,(esymm%nsym_trgk)) 387 ABI_MALLOC(esymm%tr_g0,(3,esymm%nsym_trgk)) 388 389 ! Important NOTE: 390 ! If nonsymmorphic_at_zoneborder symmetry analysis cannot be performed unless 391 ! an external database retrieved from the bilbao server (REPRES) is found. 392 idx_gk=0; idx_trgk=0 393 esymm%sgk2symrec=-999; esymm%tr_sgk2symrec=-999 394 do isym=1,nsym_fm 395 396 if (ksym_table(4,1,isym)==1) then ! S k = k +G0 397 idx_gk=idx_gk+1 398 sgk(:,:,idx_gk)=symrec_fm(:,:,isym) 399 esymm%g0(:,idx_gk)=ksym_table(1:3,1,isym) 400 esymm%sgk2symrec(idx_gk)=fm2symrec(isym) 401 if (ANY(ksym_table(1:3,1,isym)/=0).and.(ANY(ABS(Cryst%tnons(:,fm2symrec(isym)))>tol6))) then 402 esymm%nonsymmorphic_at_zoneborder=.TRUE. 403 end if 404 end if 405 406 if (esymm%can_use_tr.and.ksym_table(4,2,isym)==1) then ! -S k = k +G0 407 idx_trgk=idx_trgk+1 408 tr_sgk(:,:,idx_trgk)=symrec_fm(:,:,isym) 409 esymm%tr_g0(:,idx_trgk)=ksym_table(1:3,2,isym) 410 esymm%tr_sgk2symrec(idx_trgk)=fm2symrec(isym) 411 end if 412 end do 413 414 ABI_FREE(ksym_table) 415 ABI_FREE(symrec_fm) 416 ABI_FREE(fm2symrec) 417 418 ! ========================================== 419 ! ==== Divide the operations in classes ==== 420 ! ========================================== 421 ABI_MALLOC(dum_symafm,(esymm%nsym_gk)) 422 dum_symafm=1 423 424 !Check group closure 425 call sg_multable(esymm%nsym_gk,dum_symafm,sgk,grp_ierr) 426 ABI_CHECK(grp_ierr==0,"sg_multable failed") 427 ABI_FREE(dum_symafm) 428 429 ABI_MALLOC(tmp_nelements,(esymm%nsym_gk)) 430 ABI_MALLOC(elements_idx,(esymm%nsym_gk,esymm%nsym_gk)) 431 432 call get_classes(esymm%nsym_gk,sgk,esymm%nclass,tmp_nelements,elements_idx) 433 434 ABI_MALLOC(esymm%nelements,(esymm%nclass)) 435 esymm%nelements = tmp_nelements(1:esymm%nclass) 436 ABI_FREE(tmp_nelements) 437 438 ! From the list of symmetry operations and the lattice vectors, determine the 439 ! Bravais information including the holohedry, the centering, the coordinate of 440 ! the primitive vectors in the conventional vectors, as well as the point group, 441 msym=192; if (allocated(Cryst%symrec)) msym=size(Cryst%symrec,3) 442 ABI_MALLOC(tmp_symrec,(3,3,msym)) 443 tmp_symrec(:,:,1:esymm%nsym_gk)=sgk 444 445 call symbrav(bravais,msym,esymm%nsym_gk,ptgroup,Cryst%gprimd,tmp_symrec,tolsym,axis=sym_axis) 446 447 ABI_FREE(tmp_symrec) 448 449 write(std_out,'(a)')" symptgroup returned point group: "//TRIM(ptgroup) 450 write(std_out,'(a,i2)')" iholohedry ",bravais(1) 451 write(std_out,'(a,i2)')" center ",bravais(2) 452 write(std_out,'(a,9i3)')" gprimd in the axes of the conventional bravais lattice (*2 if center/=0)",bravais(3:11) 453 write(std_out,'(a,3i3)')" sym_axis ",sym_axis 454 455 ! Branching: 456 ! 1) If the little group is not symmorphic_at_zoneborder we can 457 ! classify the states using the irreducible representation of the point group. 458 ! 459 ! 2) If the little group is symmorphic_at_zoneborder, we have to rely on 460 ! an external database retrieved from the Bilbao server in order to classify the states. 461 ! If the file is not available, we only know the number of classes but neither their 462 ! character nor the dimension of the irreducible representation. 463 ! 464 if (esymm%nonsymmorphic_at_zoneborder) then 465 466 spgroup=0 467 chkprim=1 ! Cell must be primitive. 468 !call symlatt(bravais,std_out,msym,nptsym,ptsymrel,rprimd,tolsym) 469 !call symspgr(bravais,Cryst%nsym,spgroup,Cryst%symrel,Cryst%tnons,tolsym) 470 471 !call symanal(bravais,chkprim,genafm,msym,nsym,ptgroupma,rprimd,spgroup,symafm,symrel,tnons,tolsym) 472 473 call int2char10(spgroup,spgroup_str) 474 lgroup_fname = "lgroup_"//TRIM(spgroup_str) 475 476 if (file_exists(lgroup_fname)) then 477 ABI_ERROR("Not coded") 478 479 ! Read little groups from the external database. 480 !% call init_groupk_from_file(Lgrp,spgroup,lgroup_fname,ierr) 481 482 ! Save the irreducible representations in esymm. 483 ! Reorder symmetries such that they correspond to the Bilbao database. 484 !% allocate(esymm%Ref_irreps(esymm%nclass)) 485 !% call copy_irrep(Irreps, esymm%Ref_irreps) 486 487 else 488 write(msg,'(7a)')& 489 & "Non-symmorphic small group and zone border. ",ch10,& 490 & "External file: ",TRIM(lgroup_fname)," containing Bilbao tables not found ",ch10,& 491 & "Character analysis cannot be performed. Accidental degeneracies cannot be detected. " 492 ABI_WARNING(msg) 493 494 esymm%has_chtabs = .FALSE. 495 496 ! Reorder indices such that symmetries are packed in classes. 497 ABI_MALLOC(new_idx,(esymm%nsym_gk)) 498 ABI_MALLOC(new_g0,(3,esymm%nsym_gk)) 499 new_g0=0; iord = 0 500 do icls=1,esymm%nclass 501 do iel=1,esymm%nelements(icls) 502 iord = iord+1 503 jsym = elements_idx(iel,icls) 504 new_idx(iord) = esymm%sgk2symrec(jsym) 505 new_g0(:,iord) = esymm%g0(:,jsym) 506 end do 507 end do 508 509 esymm%sgk2symrec = new_idx 510 esymm%g0 = new_g0 511 512 ABI_FREE(new_idx) 513 ABI_FREE(new_g0) 514 end if ! file exists 515 516 else 517 ! 518 ! **** This part is still under development. It might not work for particular **** 519 ! **** orientations of the unit cell or particular lattices. **** 520 ! 521 ! The symmetries in the Bilbao database refer to the conventional unit cells. 522 ! Therefore we have to map the abinit symmetries (in reduced coordinates) 523 ! onto the Bilbao dataset. Bilbao standard settings are: 524 ! 525 ! * unique axis b (cell choice 1) for space groups withing the monoclinic system 526 ! * obverse triple hexagonal unit cell R space groups. 527 ! * origin choice two - inversion center at (0, 0, 0) - for the centrosymmetric 528 ! space groups for which there are two origins choices, within the 529 ! orthorombic, tetragonal and cubic system. 530 531 ! 1) Retrieve the rotation matrices and the irreducible representations (Bilbao setting). 532 call point_group_init(Ptg,ptgroup) 533 534 esymm%has_chtabs = .TRUE. 535 ABI_CHECK(esymm%nclass==Ptg%nclass,"esymm%nclass/=Ptg%nclass!") 536 537 do icls=1,esymm%nclass ! FIXME this is awful, should be done in a cleaner way. 538 esymm%nelements(icls)=Ptg%class_ids(2,icls) - Ptg%class_ids(1,icls) + 1 539 end do 540 541 ! 2) Generate the symmetry operations in the conventional vector coordinates. 542 conv_gprimd(:,1)=bravais(3:5) 543 conv_gprimd(:,2)=bravais(6:8) 544 conv_gprimd(:,3)=bravais(9:11) 545 546 axes = conv_gprimd 547 call matr3inv(conv_gprimd,axes) !; axes=TRANSPOSE(axes) 548 549 conv_gprimd=MATMUL(Cryst%gprimd,TRANSPOSE(axes)) 550 !conv_gprimd=MATMUL(axes,Cryst%gprimd) 551 !conv_gprimd=MATMUL(TRANSPOSE(axes),Cryst%gprimd) 552 !write(std_out,*)"conv_gprimd:", conv_gprimd 553 554 ptgroup_name = ADJUSTL(ptgroup) 555 556 select case (ptgroup_name) 557 558 case ("3m","-3m") 559 call wrtout(std_out," Changing the conventional cell: rhombohedral --> triple hexagonal","COLL") 560 ! Transformation matrices: primitive rhombohedral --> triple hexagonal cell obverse setting. Table 5.1.3.1 ITA page 81. 561 pmat1 = RESHAPE( (/ 1,-1, 0, 0, 1,-1, 1, 1, 1/), (/3,3/) ) ! R1 562 pmat2 = RESHAPE( (/ 0, 1,-1,-1, 0, 1, 1, 1, 1/), (/3,3/) ) ! R2 563 pmat3 = RESHAPE( (/-1, 0, 1, 1,-1, 0, 1, 1, 1/), (/3,3/) ) ! R3 564 pmat4 = RESHAPE( (/-1, 1, 0, 0,-1, 1, 1, 1, 1/), (/3,3/) ) ! R1 reverse setting. 565 pmat5 = RESHAPE( (/ 0,-1, 1, 1, 0,-1, 1, 1, 1/), (/3,3/) ) ! R2 reverse setting. 566 pmat6 = RESHAPE( (/ 1, 0,-1,-1, 1, 0, 1, 1, 1/), (/3,3/) ) ! R3 reverse setting. 567 conv_gprimd = MATMUL(conv_gprimd,pmat1) 568 !conv_gprimd = MATMUL(conv_gprimd,pmat2) 569 !conv_gprimd = MATMUL(conv_gprimd,pmat3) 570 !conv_gprimd = MATMUL(conv_gprimd,pmat4) 571 !conv_gprimd = MATMUL(conv_gprimd,pmat5) 572 !conv_gprimd = MATMUL(conv_gprimd,pmat6) 573 !write(std_out,*)" New conv_gprimd:", conv_gprimd 574 575 case ("mm2") 576 call wrtout(std_out," Changing the conventional cell: unconventional orthorhombic setting --> conventional","COLL") 577 ! Transformation matrices: unconvential orthorhombic --> conventional orthorhombic. Table 5.1.3.1 ITA page 81. 578 pmat1 = RESHAPE( (/ 0, 1, 0, 1, 0, 0, 0, 0,-1/), (/3,3/) ) ! ( b, a,-c) --> (a,b,c) 579 pmat2 = RESHAPE( (/ 0, 1, 0, 0, 0, 1, 1, 0, 0/), (/3,3/) ) ! ( c, a, b) --> (a,b,c) 580 pmat3 = RESHAPE( (/ 0, 0, 1, 0, 1, 0,-1, 0, 0/), (/3,3/) ) ! (-c, b, a) --> (a,b,c) 581 pmat4 = RESHAPE( (/ 0, 0, 1, 1, 0, 0, 0, 1, 0/), (/3,3/) ) ! ( b, c, a) --> (a,b,c) 582 pmat5 = RESHAPE( (/ 1, 0, 0, 0, 0, 1, 0,-1, 0/), (/3,3/) ) ! ( a,-c, b) --> (a,b,c) 583 conv_gprimd = MATMUL(conv_gprimd,pmat2) 584 !write(std_out,*)" New conv_gprimd:", conv_gprimd 585 case default 586 continue 587 end select 588 589 ABI_MALLOC(conv_symrec,(3,3,esymm%nsym_gk)) 590 conv_symrec = sgk 591 592 !axes=zero; axes(1,1)=one ; axes(2,2)=one ; axes(3,3)=one 593 !call symrelrot(esymm%nsym_gk,conv_gprimd,axes,conv_symrec,tolsym) 594 call symrelrot(esymm%nsym_gk,Cryst%gprimd,conv_gprimd,conv_symrec,tolsym) 595 596 ! 3) Reorder indices such that symmetries are packed in classes. 597 ABI_MALLOC(found,(esymm%nsym_gk)) 598 ABI_MALLOC(new_idx,(esymm%nsym_gk)) 599 ABI_MALLOC(new_g0,(3,esymm%nsym_gk)) 600 new_g0=0; found=0 601 602 do isym=1,esymm%nsym_gk 603 do jsym=1,esymm%nsym_gk 604 if (ALL(Ptg%sym(:,:,isym) == conv_symrec(:,:,jsym) )) then 605 found(isym) = found(isym) + 1 606 new_idx(isym) = esymm%sgk2symrec(jsym) 607 new_g0(:,isym) = esymm%g0(:,jsym) 608 !EXIT 609 end if 610 end do 611 end do 612 ! 613 ! DEBUGGING SECTION 614 !do isym=1,esymm%nsym_gk 615 ! jsym=esymm%sgk2symrec(isym) 616 ! call print_symmetries(1,Cryst%symrec(:,:,jsym),Cryst%tnons(:,jsym),Cryst%symafm(jsym)) 617 ! write(std_out,*)esymm%g0(:,isym) 618 !end do 619 620 if ( Ptg%nsym/=esymm%nsym_gk .or. ANY(found/=1) ) then 621 !write(std_out,*)Ptg%nsym, esymm%nsym_gk 622 !write(std_out,'(a,(i2))')" found = ",found 623 write(std_out,*)" Ptg%sym list, conv_symrec list, found Ptg% " 624 do isym=1,Ptg%nsym 625 write(std_out,'(a,i2,a,9i2,4x,a,9i2)')" found ",found(isym)," Ptg ",Ptg%sym(:,:,isym),"conv_symrec ",conv_symrec(:,:,isym) 626 end do 627 msg = " sgk and esymm%Ptg are inconsistent. Check tables or source" 628 ABI_WARNING(msg) 629 esymm%err_msg = msg(1:500) 630 esymm%err_status = ESYM_PTG_WRONG_MAPPING 631 esymm%has_chtabs = .FALSE. 632 633 else ! Reorder symmetries. 634 esymm%sgk2symrec = new_idx 635 esymm%g0 = new_g0 636 end if 637 638 ABI_FREE(new_idx) 639 ABI_FREE(new_g0) 640 ABI_FREE(found) 641 ABI_FREE(conv_symrec) 642 643 if (esymm%has_chtabs) then 644 ! Multiply the point group irreps by e^{-ik.\tau} to have the irreps of the little group. 645 ! Store the results in esymm%Ref_irreps so that one can classify the states afterwards. 646 ABI_MALLOC(esymm%Ref_irreps,(esymm%nclass)) 647 ABI_MALLOC(phase_mkt,(esymm%nsym_gk)) 648 649 do isym_gk=1,esymm%nsym_gk 650 isym = esymm%sgk2symrec(isym_gk) 651 mkt = -two_pi * DOT_PRODUCT(esymm%kpt, Cryst%tnons(:,isym)) 652 phase_mkt(isym_gk) = CMPLX(DCOS(mkt), DSIN(mkt)) 653 end do 654 655 call copy_irrep(Ptg%Irreps,esymm%Ref_irreps,phase_mkt) 656 ABI_FREE(phase_mkt) 657 end if 658 659 #if 0 660 ! Herring test requires the evaluation of the expression: 661 ! 662 ! sum_{S,\tau} \chi^{k,\alpha} ({S|\tau}^2) 663 ! 664 ! where Sk = -k + g0, and \chi is the trace of the \alpha-th 665 ! irreducible representation of the little group of k. 666 ! \chi^{k,\alpha} = e^{-ik.\tau} \chi(\alpha) provided that 667 ! we are not at zone border with a non-symmorphic operation. 668 ! The expression is always real and it can only be equal to \pm Ptg%nsym or zero. 669 ! FIXME this part has to be rewritten from scratch. 670 !if (esymm%err_status/=esymm_NOERROR) then 671 ! write(std_out,*)" Skipping Herring test" 672 ! goto 110 673 !end if 674 675 if (esymm%can_use_tr) then 676 ABI_MALLOC(her_test,(esymm%nclass)) 677 678 ABI_MALLOC(tr_conv_symrec,(3,3,esymm%nsym_trgk)) 679 do isym_tr=1,esymm%nsym_trgk 680 isym = esymm%tr_sgk2symrec(isym_tr) 681 tr_conv_symrec(:,:,isym_tr)=Cryst%symrec(:,:,isym) 682 end do 683 684 call symrelrot(esymm%nsym_trgk,Cryst%gprimd,conv_gprimd,tr_conv_symrec_tr,tolsym) 685 686 do isym_tr=1,esymm%nsym_trgk 687 isym = esymm%tr_sgk2symrec(isym_tr) 688 !rot2 = MATMUL(tr_sgk(:,:,isym),tr_sgk(:,:,isym)) 689 !tau2 = MATMUL(tr_sgk(:,:,isym),Cryst%tnons(:,isym)) + Cryst%tnons(:,isym) 690 691 rot2 = MATMUL(tr_conv_symrec(:,:,isym_tr),tr_conv_symrec(:,:,isym_tr)) 692 tau2 = MATMUL(tr_conv_symrec(:,:,isym_tr),Cryst%tnons(:,isym)) + Cryst%tnons(:,isym) 693 694 phase_k = EXP(-j_dpc*two_pi*DOT_PRODUCT(kpoint,tau2)) 695 call locate_sym(Ptg,rot2,isym2,icls2) 696 697 do irp=1,esymm%nclass 698 her_test(irp) = her_test(irp) + phase_k * Ptg%Irreps(irp)%trace(icls2) 699 end do 700 end do 701 702 ABI_FREE(tr_conv_symrec) 703 704 ! FIXME 705 ABI_MALLOC(esymm%herring_test,(esymm%nclass)) 706 707 do irp=1,esymm%nclass 708 if ( ABS(her_test(irp) - Ptg%nsym) < tol6 ) then 709 esymm%herring_test(irp) = +1 710 else if ( ABS(her_test(irp)) < tol6 ) then 711 esymm%herring_test(irp) = 0 712 else if ( ABS(her_test(irp) + Ptg%nsym) < tol6 ) then 713 esymm%herring_test(irp) = -1 714 else 715 write(msg,'(a,i2,2a,i0,a,i2)')& 716 & "Herring test for the irreducible representation number ",irp,ch10,& 717 & "gave ",esymm%herring_test(irp),", while it should be 0 or +- ",Ptg%nsym 718 ABI_WARNING(msg) 719 esymm%err_msg =msg 720 esymm%err_status=esymm_HERRING_WRONG_TEST 721 end if 722 end do 723 724 ABI_FREE(her_test) 725 end if ! can_use_tr 726 #endif 727 ! 728 ! Final check 729 !allocate(mtab(esymm%nsym_gk,esymm%nsym_gk)) 730 !call mult_table(esymm%nsym_gk,Ptg%sym,mtab) 731 732 !do isym=1,esymm%nsym_gk 733 ! isym1 = esymm%sgk2symrec(isym) 734 ! do jsym=1,esymm%nsym_gk 735 ! isym2 = esymm%sgk2symrec(jsym) 736 ! rot2 = MATMUL(Cryst%symrec(:,:,isym1),Cryst%symrec(:,:,isym2)) 737 738 ! iprod = mtab(isym,jsym) 739 740 ! do irp=1,esymm%nclass 741 ! dim_irrep = Ptg%Irreps(irp)%dim 742 ! allocate(mat_test(dim_irrep,dim_irrep)) 743 ! mat_test = Ptg%Irreps(irp)%mat(:,:,isym) * Ptg%Irreps(irp)%mat(:,:,jsym) 744 ! !call locate_sym(Ptg,rot2,isym2,icls2) 745 ! write(std_out,*)mat_test - Ptg%Irreps(irp)%mat(:,:,iprod) 746 ! deallocate(mat_test) 747 ! end do 748 ! 749 ! end do 750 !end do 751 ! 752 !deallocate(mtab) 753 end if 754 755 ABI_FREE(sgk) 756 ABI_FREE(tr_sgk) 757 ABI_FREE(elements_idx) 758 759 !% allocate(esymm%irrep2b(0:esymm%nclass)) 760 !% call nullify_coeff(esymm%irrep2b) 761 ! 762 ! 1) Allocate space for the irreducible representations. 763 764 ! 2) Try to determine if we are in presence of an accidental degeneracy. Sufficient condition: 765 ! There exists a set of degenerate states whose dimension is greater than the dimension 766 ! of the irreducible representations of the point group. The check can be done only 767 ! if Character tables are available. 768 769 if (esymm%has_chtabs) then 770 ABI_MALLOC(dim_irreps,(esymm%nclass)) 771 dim_irreps = (/(esymm%Ref_irreps(irp)%dim, irp=1,esymm%nclass)/) 772 end if 773 774 nacc_deg=0 775 ABI_MALLOC(esymm%degs_dim,(esymm%ndegs)) 776 ABI_MALLOC(esymm%Calc_irreps,(esymm%ndegs)) 777 778 if (esymm%can_use_tr) then 779 ABI_MALLOC(esymm%trCalc_irreps,(esymm%ndegs)) 780 end if 781 782 do idg=1,esymm%ndegs 783 dim_degs=esymm%degs_bounds(2,idg)-esymm%degs_bounds(1,idg)+1 784 785 if (esymm%has_chtabs) then 786 if (ALL(dim_degs /= dim_irreps)) then ! An accidental degeneracy is present. 787 nacc_deg=nacc_deg+1 788 end if 789 end if 790 791 esymm%degs_dim(idg) = dim_degs 792 793 call init_irrep(esymm%Calc_irreps(idg),esymm%nsym_gk,dim_degs) 794 if (esymm%can_use_tr) then 795 call init_irrep(esymm%trCalc_irreps(idg),esymm%nsym_trgk,dim_degs) 796 end if 797 end do ! idg 798 799 if (esymm%has_chtabs) then 800 ABI_FREE(dim_irreps) 801 if (nacc_deg/=0) then 802 write(msg,'(a,i0,a)')" Detected ",nacc_deg," accidental degeneracies." 803 ABI_WARNING(msg) 804 esymm%err_status=ESYM_ACCDEG_ERROR 805 ! TODO this should signal to the caller that we have to decompose the calculated representation. 806 esymm%err_msg =msg(1:500) 807 end if 808 end if 809 810 DBG_EXIT("COLL") 811 812 end subroutine esymm_init
m_esymm/esymm_print [ Functions ]
[ Top ] [ m_esymm ] [ Functions ]
NAME
esymm_print
FUNCTION
INPUTS
OUTPUT
SOURCE
829 subroutine esymm_print(esymm,unit,mode_paral,prtvol) 830 831 !Arguments ------------------------------------ 832 !scalars 833 integer,optional,intent(in) :: prtvol,unit 834 character(len=4),optional,intent(in) :: mode_paral 835 type(esymm_t),intent(in) :: esymm 836 837 !Local variables------------------------------- 838 !scalars 839 integer :: icl,idg,my_unt,my_prtvol 840 integer :: irr_idx,nstates,nunknown,istart,istop,ii 841 character(len=4) :: my_mode 842 character(len=1000) :: fmt,msg,msg0 843 844 ! ********************************************************************* 845 846 my_unt =std_out; if (PRESENT(unit )) my_unt =unit 847 my_prtvol=0 ; if (PRESENT(prtvol )) my_prtvol=prtvol 848 my_mode ='COLL' ; if (PRESENT(mode_paral)) my_mode =mode_paral 849 850 write(fmt,*)'(2a,3f8.4,3a,i4,2a,i3,2a,i2,2a,i2,a,',esymm%nclass,'i2,a)' 851 write(msg,fmt)ch10,& 852 & ' ===== Character of bands at k-point: ',esymm%kpt,' ===== ',ch10,& 853 & ' Total number of bands analyzed .................. ',esymm%nbnds,ch10,& 854 & ' Number of degenerate sets detected .............. ',esymm%ndegs,ch10,& 855 & ' Number of operations in the little group of k ... ',esymm%nsym_gk,ch10,& 856 & ' Number of classes (irreps) in the group of k .... ',esymm%nclass,' (',(esymm%nelements(icl),icl=1,esymm%nclass),' )' 857 call wrtout(my_unt,msg,my_mode) 858 859 if (esymm%nonsymmorphic_at_zoneborder) then 860 call wrtout(my_unt," Non-symmorphic small group at zone border. Character analysis not available ",my_mode) 861 end if 862 863 if (esymm_failed(esymm)) then 864 write(std_out,'(3a)')"Band classification algorithm failed with the error:",ch10,TRIM(esymm%err_msg) 865 write(msg,'(3a)')"Band classification algorithm failed with the error:",ch10,TRIM(esymm%err_msg) 866 call wrtout(my_unt,msg,my_mode) 867 end if 868 869 !nunknown=0 870 !do iband=1,esymm%nbnds 871 ! irr_idx = esymm%b2irrep(iband) 872 ! if (irr_idx /= 0) then 873 ! if ( esymm%has_chtabs) irr_name = esymm%Ref_Irreps(irr_idx)%name 874 ! if (.not.esymm%has_chtabs) write(irr_name,'(i0)')irr_idx ! use the index instead of the name. 875 ! else 876 ! irr_name = "???" 877 ! nunknown = nunknown +1 878 ! end if 879 ! write(msg,'(a,i3,2a)')' Band ',iband,' belongs to irrep ',TRIM(irr_name) 880 ! call wrtout(my_unt,msg,my_mode) 881 !end do 882 883 do irr_idx=1,esymm%nclass 884 nstates = size(esymm%irrep2b(irr_idx)%value) 885 if (esymm%has_chtabs) then 886 write(msg0,'(a,i0,3a)')" Found ",nstates," states with character ",TRIM(esymm%Ref_irreps(irr_idx)%name),": " 887 else 888 write(msg0,'(2(a,i0),a)')" Found ",nstates," states with character index ",irr_idx,": " 889 end if 890 do istart=1,nstates,20 891 istop=istart+11; if (istop>nstates) istop=nstates 892 write(msg,'(20(1x,i0))')(esymm%irrep2b(irr_idx)%value(ii), ii=istart,istop) 893 if (istart==1) msg = TRIM(msg0)//TRIM(msg) 894 if (istart/=1) msg = " "//TRIM(msg) 895 call wrtout(my_unt,msg,my_mode) 896 end do 897 end do 898 899 nunknown = size(esymm%irrep2b(0)%value) 900 if (nunknown > 0) then 901 write(msg0,'(a,i0,a)')" WARNING: ",nunknown," states have not been classified:" 902 do istart=1,nunknown,20 903 istop=istart+11; if (istop>nunknown) istop=nunknown 904 write(msg,'(20(1x,i0))')(esymm%irrep2b(0)%value(ii), ii=istart,istop) 905 if (istart==1) msg = TRIM(msg0)//TRIM(msg) 906 if (istart/=1) msg = " "//TRIM(msg) 907 call wrtout(my_unt,msg,my_mode) 908 end do 909 end if 910 911 if (my_prtvol>0 .or. nunknown>0 .or. .not.esymm%has_chtabs) then ! print the calculated character table. 912 call wrtout(my_unt,ch10//" Calculated character table ",my_mode) 913 !write(fmt,*)'(i2,a,i2,1x,',esymm%nclass,'(a,2f6.3),a)' 914 write(fmt,*)'(i2,a,i2,1x,',esymm%nclass,'(a,2f5.2),a)' 915 do idg=1,esymm%ndegs 916 write(msg,fmt)& 917 & esymm%degs_bounds(1,idg),'-',esymm%degs_bounds(2,idg),& 918 & ('|',esymm%Calc_irreps(idg)%trace(esymm%nelements(icl)), icl=1,esymm%nclass),'|' 919 call wrtout(my_unt,msg,my_mode) 920 end do 921 end if 922 923 end subroutine esymm_print
m_esymm/esymm_symmetrize_mels [ Functions ]
[ Top ] [ m_esymm ] [ Functions ]
NAME
esymm_symmetrize_mels
FUNCTION
INPUTS
esymm<esymm_t>
SOURCE
1311 subroutine esymm_symmetrize_mels(esymm,lbnd,ubnd,in_me,out_me) 1312 1313 !Arguments ------------------------------------ 1314 !scalars 1315 integer :: lbnd,ubnd 1316 type(esymm_t),target,intent(in) :: esymm 1317 !arrays 1318 complex(dpc),intent(in) :: in_me(2,lbnd:ubnd,lbnd:ubnd) 1319 complex(dpc),intent(out) :: out_me(lbnd:ubnd,lbnd:ubnd) 1320 1321 !Local variables------------------------------- 1322 !scalars 1323 integer :: idg1,b1_start,b1_stop,irp1 1324 integer :: idg2,b2_start,b2_stop,irp2 1325 integer :: ii,jj,ib,jb,kk,kb,lb,ll 1326 complex(dpc) :: tr_ofd,ofd,dsd,tr_dsd 1327 type(irrep_t),pointer :: Irrep1,Irrep2 1328 type(irrep_t),pointer :: tr_Irrep1,tr_Irrep2 1329 1330 ! ********************************************************************* 1331 1332 if (esymm_failed(esymm)) then 1333 ABI_ERROR("Symmetrization cannot be performed. You should not be here!") 1334 end if 1335 1336 do idg1=1,esymm%ndegs ! First loop over set of degenerate states. 1337 b1_start = esymm%degs_bounds(1,idg1) 1338 b1_stop = esymm%degs_bounds(2,idg1) 1339 1340 !if (b1_stop<lbnd .or. b2_start >ubnd) then 1341 ! ABI_ERROR("Wrong band indices, check esymm initialization") 1342 !end if 1343 1344 Irrep1 => esymm%Calc_irreps(idg1) 1345 if (esymm%can_use_tr) tr_Irrep1 => esymm%trCalc_irreps(idg1) 1346 irp1 = esymm%b2irrep(b1_start) 1347 1348 do idg2=1,esymm%ndegs ! Second loop over set of degenerate states. 1349 !write(std_out,*)" ==> Symmetrizing degenerate set ",idg1,idg2 1350 b2_start = esymm%degs_bounds(1,idg2) 1351 b2_stop = esymm%degs_bounds(2,idg2) 1352 irp2 = esymm%b2irrep(b2_start) 1353 1354 if (irp1/=irp2 .or. idg1==idg2) CYCLE ! Skip diago elements or elements belonging to different irreps. 1355 1356 Irrep2 => esymm%Calc_irreps(idg2) 1357 if (esymm%can_use_tr) tr_Irrep2 => esymm%trCalc_irreps(idg2) 1358 ! 1359 ! Symmetrize the off-diagonal matrix elements. 1360 ! summing over kk and ll. ii and jj are the indices of the bands that are symmetrized 1361 do ii=1,b1_stop-b1_start+1 1362 ib= ii+b1_start-1 1363 do jj=1,b2_stop-b2_start+1 1364 jb= jj+b2_start-1 1365 !write(std_out,*)" ====> Symmetrizing ",ib,jb 1366 1367 ofd= czero; tr_ofd=czero 1368 do kk=1,b1_stop-b1_start+1 1369 kb= kk+b1_start-1 1370 do ll=1,b2_stop-b2_start+1 1371 lb= ll+b2_start-1 1372 dsd = sum_irreps(Irrep1,Irrep2,kk,ii,ll,jj) 1373 ofd = ofd + dsd * in_me(1,kb,lb) 1374 if (esymm%can_use_tr) then 1375 tr_dsd = sum_irreps(tr_Irrep1,tr_Irrep2,kk,jj,ll,ii) ! Exchange of band indices. 1376 tr_ofd = tr_ofd + tr_dsd * in_me(2,kb,lb) ! Contribution obtained from TR. 1377 end if 1378 end do 1379 end do 1380 1381 out_me(ib,jb)= ofd/esymm%nsym_gk 1382 if (esymm%can_use_tr .and. esymm%nsym_trgk>0) out_me(ib,jb)= out_me(ib,jb) + tr_ofd/esymm%nsym_trgk 1383 end do 1384 end do 1385 end do 1386 end do 1387 1388 end subroutine esymm_symmetrize_mels
m_esymm/esymm_t [ 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
68 type,public :: esymm_t 69 70 integer :: nspinor 71 ! Number of spinorial components. 72 73 integer :: first_ib 74 ! Index of the first treated band. 75 76 integer :: nbnds 77 ! Number of bands for this k-point and spin. 78 79 integer :: nclass 80 ! The number of classes in the group of k. 81 82 integer :: nsym_gk 83 ! Number of symmetries in the group of k. Namely that the set of symmetries such that Sk = k +G0. 84 85 integer :: nsym_trgk 86 ! Number of symmetries in the extended group of k. Namely that the set of symmetries such that -Sk = k + G0. 87 88 integer :: err_status = ESYM_NOERROR 89 ! Flag signaling if the classification algorithm succeed or not. 90 91 real(dp) :: tol_deg 92 ! Energy tolerance below which two states are considered degenerate. 93 94 logical :: can_use_tr 95 ! .TRUE. if time-reversal can be used 96 97 logical :: only_trace 98 ! if .TRUE. only the trace of a single matrix per class is calculated 99 ! this is the standard way used to analyze bands symmetries. If .FALSE. 100 ! the full matrices of the irreducible representations are calculated and stored 101 102 logical :: has_spatial_inv 103 ! .TRUE. if the inversion belongs to the space group 104 105 logical :: nonsymmorphic_at_zoneborder 106 ! if .TRUE. analysis cannot be performed since kpt is 107 ! at border zone and non-zero fractional translations are present in the space group 108 109 logical :: has_chtabs 110 ! True if Ref_irreps and character tables are available (tables are initialized either 111 ! from point group irreps or from an external database downloaded from the Bilbao server) 112 113 real(dp) :: kpt(3) 114 ! The crystalline momentum of the wavefunctions in reduced coordinates. 115 116 character(len=500) :: err_msg="None" 117 118 integer,allocatable :: g0(:,:) 119 ! g0(3,nsym_gk) 120 ! The umklapp g0 vector associated to each little group operation. 121 122 integer,allocatable :: tr_g0(:,:) 123 ! tr_g0(3,nsym_trgk) 124 ! The umklapp g0 vector associated to each little group operation. 125 126 integer :: ndegs 127 ! Number of degenerate states. 128 129 integer,allocatable :: nelements(:) 130 ! nelements(nclass) 131 ! Number of symmetry operations in each class. 132 133 integer,allocatable :: sgk2symrec(:) 134 ! sgk2symrec(nsym_gk) 135 ! Mapping between the symmetries of the group of k and the symrec(l) array. 136 ! The symmetries of the little group are always packed in classes to facilitate 137 ! the calculation of the character of the irrep. Abinit symmetries are randomly ordered. 138 139 integer,allocatable :: tr_sgk2symrec(:) 140 ! trsgk2symrec(nsym_trgk) 141 ! Mapping between the symmetries of the group of k and the symrec(l) array. 142 ! The symmetries of the little group are always packed in classes to facilitate 143 ! the calculation of the character of the irrep. Abinit symmetries are randomly ordered. 144 145 integer,allocatable :: herring_test(:) 146 ! herring_test(nclass) 147 ! The result of Herring test for each irreducible representantion of the group of k. 148 ! Possible values are: +1, 0, -1 149 150 integer,allocatable :: b2irrep(:) 151 ! b2irrep(nbnds) 152 ! For each band, it gives the index of the irreducible representation in Ref_irreps. 153 154 type(coeffi1_type),allocatable :: irrep2b(:) 155 ! irrep2b(0:nclass)%value(:) 156 ! Ragged arrays with the mapping between the set of irreducible representation and the band indices. 157 ! irrep2b(irp)%value(:) gives the indices of the states belonging to irrep irp, irp=1,nclass 158 ! irrep2b(0)%value(:) stores the indices of the states that have not been classified due to 159 ! the presence of an accidental degeneracy. 160 161 integer,allocatable :: degs_bounds(:,:) 162 ! degs_bounds(2,ndegs) 163 ! degs_bounds(1,idg)= first band index of the degenerate set idg=1,ndegs 164 ! degs_bounds(2,idg)= final band index of the degenerate set idg=1,ndegs 165 166 integer,allocatable :: degs_dim(:) 167 ! degs_dim(ndegs) 168 ! Number of states in each degenerate subspace. Cannot be larger that nclass provided 169 ! that no accidental degeneracy occurs. 170 171 !% integer,allocatable :: class_ids(:,:) 172 ! class_ids(2,nclass) 173 ! (1,icl) = index of the first symmetry of class icl 174 ! (2,icl) = index of the last symmetry of class icl 175 ! Note that symmetries in sym are packed in classes. 176 177 type(irrep_t),allocatable :: Calc_irreps(:) 178 ! Calc_irreps(ndegs) 179 ! The representations of the little group of k calculated from the wavefunctions. <\phi_nk|R_t|\phi_mk> 180 ! where R_t belong to the little group of k. 181 ! They represent an unitary irreducible representation provided that no accidental degeneracy occurs. 182 183 type(irrep_t),allocatable :: trCalc_irreps(:) 184 ! trCalc_irreps(ndegs) 185 ! The representations of the little group of k calculated from the wavefunctions. <\phi_nk|R_t|\phi_mk> 186 ! where R_t belong to the little group of k. 187 ! They represent an unitary irreducible representation provided that no accidental degeneracy occurs. 188 189 type(irrep_t),allocatable :: Ref_irreps(:) 190 ! Irreps(nclass) 191 ! Reference irreducible representations of the group of k derived from the point group 192 ! or from the external database downloaded from the Bilbao web site. 193 194 end type esymm_t 195 196 public :: esymm_init ! Initialize the object 197 public :: esymm_print ! Print info 198 public :: esymm_free ! Free memory 199 public :: esymm_finalize ! Finalize the object 200 public :: esymm_symmetrize_mels ! Symmetrize given matrix elements 201 public :: esymm_failed ! True if symmetry analysis failed.
m_esymm/polish_irreps [ Functions ]
[ Top ] [ m_esymm ] [ Functions ]
NAME
polish_irreps
FUNCTION
INPUTS
SOURCE
1430 subroutine polish_irreps(Irreps) 1431 1432 !Arguments ------------------------------------ 1433 !scalars 1434 type(irrep_t),intent(inout) :: Irreps(:) 1435 1436 !Local variables------------------------------- 1437 !scalars 1438 integer,parameter :: ldvl1=1,ldvr1=1 1439 integer :: irp,sym,dim,ldvr,ii,ivec,jvec,info 1440 !character(len=500) :: msg 1441 !arrays 1442 complex(dpc),allocatable :: vl(:,:),vr(:,:),vrm1(:,:),overlap(:,:) 1443 complex(dpc),allocatable :: cmat(:,:),eigval(:) 1444 1445 ! ********************************************************************* 1446 1447 ! Eigen decomposition: A = V D V^{-1}. 1448 do irp=1,SIZE(Irreps) 1449 dim = Irreps(irp)%dim 1450 ABI_MALLOC(cmat,(dim,dim)) 1451 ABI_MALLOC(eigval,(dim)) 1452 ldvr=dim 1453 ABI_MALLOC(vl,(ldvl1,dim)) 1454 ABI_MALLOC(vr,(ldvr,dim)) 1455 ABI_MALLOC(vrm1,(dim,dim)) 1456 ABI_MALLOC(overlap,(dim,dim)) 1457 do sym=1,Irreps(irp)%nsym 1458 cmat = Irreps(irp)%mat(:,:,sym) 1459 call xgeev("No vectors","Vectors",dim,cmat,dim,eigval,vl,ldvl1,vr,ldvr) 1460 ! 1461 ! Orthogonalize the eigenvectors using Cholesky orthogonalization. 1462 do jvec=1,dim 1463 do ivec=1,jvec 1464 overlap(ivec,jvec) = DOT_PRODUCT(vr(:,ivec),vr(:,jvec)) 1465 end do 1466 end do 1467 ! 1468 ! 2) Cholesky factorization: overlap = U^H U with U upper triangle matrix. 1469 call ZPOTRF('U',dim,overlap,dim,info) 1470 ABI_CHECK(info == 0, sjoin('ZPOTRF returned info=', itoa(info))) 1471 1472 ! 3) Solve X U = Vr, on exit the Vr treated by this node is orthonormalized. 1473 call ZTRSM('R','U','N','N',dim,dim,cone,overlap,dim,vr,dim) 1474 1475 !write(std_out,*)"After ortho",MATMUL(TRANSPOSE(CONJG(vr)),vr) 1476 1477 vrm1 = vr 1478 call xginv(vrm1,dim) 1479 do ii=1,dim 1480 eigval(ii) = eigval(ii)/ABS(eigval(ii)) ! Rescale the eigevalues. 1481 vrm1(ii,:) = eigval(ii) * vrm1(ii,:) 1482 end do 1483 Irreps(irp)%mat(:,:,sym) = MATMUL(vr,vrm1) 1484 Irreps(irp)%trace(sym) = get_trace(Irreps(irp)%mat(:,:,sym)) 1485 end do 1486 ABI_FREE(cmat) 1487 ABI_FREE(eigval) 1488 ABI_FREE(vl) 1489 ABI_FREE(vr) 1490 ABI_FREE(vrm1) 1491 ABI_FREE(overlap) 1492 end do 1493 1494 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.
SOURCE
1269 function which_irrep(esymm,trace,tolerr) 1270 1271 !Arguments ------------------------------------ 1272 !scalars 1273 integer :: which_irrep 1274 real(dp),intent(in) :: tolerr 1275 type(esymm_t),intent(in) :: esymm 1276 !arrays 1277 complex(dpc),intent(in) :: trace(esymm%nsym_gk) 1278 1279 !Local variables------------------------------- 1280 !scalars 1281 integer :: irp 1282 1283 ! ********************************************************************* 1284 1285 which_irrep = 0 1286 if (esymm%has_chtabs) then ! Symmetry analysis can be performed. 1287 do irp=1,esymm%nclass 1288 if ( ALL( ABS(esymm%Ref_irreps(irp)%trace(:) - trace(:)) < tolerr)) then 1289 which_irrep = irp; EXIT 1290 end if 1291 end do 1292 end if 1293 1294 end function which_irrep