TABLE OF CONTENTS


ABINIT/complete_gamma [ Functions ]

[ Top ] [ Functions ]

NAME

 complete_gamma

FUNCTION

 Use the set of special q points calculated by the Monkhorst & Pack Technique.
 Check if all the informations for the q points are present in the input gamma matrices.
 Generate the gamma matrices (already summed over the FS) of the set of q points which
 samples homogeneously the entire Brillouin zone.

INPUTS

 qpttoqpt = qpoint index mapping under symops

OUTPUT

 gamma_qpt = in/out: set of gamma matrix elements completed and symmetrized
    gamma_qpt(2,nbranch**2,nsppol,nqpt_full)

PARENTS

      elphon,integrate_gamma_alt,m_phgamma

CHILDREN

      zgemm

SOURCE

 853 subroutine complete_gamma(Cryst,nbranch,nsppol,nqptirred,nqpt_full,ep_scalprod,qirredtofull,qpttoqpt,gamma_qpt)
 854 
 855 
 856 !This section has been created automatically by the script Abilint (TD).
 857 !Do not modify the following lines by hand.
 858 #undef ABI_FUNC
 859 #define ABI_FUNC 'complete_gamma'
 860 !End of the abilint section
 861 
 862  implicit none
 863 
 864 !Arguments ------------------------------------
 865 !scalars
 866  integer,intent(in) :: nsppol,nbranch,nqptirred,nqpt_full,ep_scalprod
 867  type(crystal_t),intent(in) :: Cryst
 868 !arrays
 869  integer,intent(in) :: qirredtofull(nqptirred)
 870  integer,intent(in) :: qpttoqpt(2,Cryst%nsym,nqpt_full)
 871  real(dp), intent(inout) :: gamma_qpt(2,nbranch**2,nsppol,nqpt_full)
 872 
 873 !Local variables-------------------------------
 874 !scalars
 875  integer :: ibranch,ieqqpt,ii,natom,nsym,iqpt,isppol,isym
 876  integer :: itim,jbranch,jj,kk,ll,neqqpt,iatom,ancestor_iatom,iqpt_fullbz
 877 !arrays
 878  integer :: symmetrized_qpt(nqpt_full)
 879  integer :: gkk_flag(nbranch,nbranch,nsppol,nqpt_full)
 880  real(dp) :: ss(3,3)
 881  real(dp) :: tmp_mat(2,nbranch,nbranch)
 882  real(dp) :: tmp_mat2(2,nbranch,nbranch)
 883  real(dp) :: ss_allatoms(2,nbranch,nbranch)
 884  real(dp) :: c_one(2), c_zero(2)
 885  real(dp),allocatable :: gkk_qpt_new(:,:,:),gkk_qpt_tmp(:,:,:)
 886 
 887 ! *********************************************************************
 888 
 889  c_one = (/one,zero/)
 890  c_zero = (/zero,zero/)
 891 
 892  natom = Cryst%natom
 893  nsym  = Cryst%nsym
 894 
 895 !Generation of the gkk matrices relative to the q points
 896 !of the set which samples the entire Brillouin zone
 897 
 898 !set up flags for gamma_qpt matrices we have
 899  gkk_flag = -1
 900  do iqpt=1,nqptirred
 901    iqpt_fullbz = qirredtofull(iqpt)
 902    gkk_flag(:,:,:,iqpt_fullbz) = 1
 903  end do
 904 
 905  symmetrized_qpt(:) = -1
 906 
 907  ABI_ALLOCATE(gkk_qpt_new,(2,nbranch**2,nsppol))
 908  ABI_ALLOCATE(gkk_qpt_tmp,(2,nbranch**2,nsppol))
 909 
 910  do iqpt=1,nqpt_full
 911 !
 912 !  Already symmetrized?
 913    if (symmetrized_qpt(iqpt) == 1) cycle
 914 
 915    gkk_qpt_new(:,:,:) = zero
 916 
 917 !  loop over qpoints equivalent to iqpt
 918    neqqpt=0
 919 !  do not use time reversal symmetry to complete the qpoints:
 920 !  do not know what happens to the gamma matrices
 921 !  11/2011: MJV: time reversal is needed here if inversion is absent
 922 !  - used in read_gkk and all reductions of q-points by symmetry.
 923 
 924    do itim=1,2
 925      do isym=1,nsym
 926 !      ieqqpt is sent onto iqpt by itim/isym
 927        ieqqpt = qpttoqpt(itim,isym,iqpt)
 928 
 929 
 930        if (gkk_flag(1,1,1,ieqqpt) == -1) cycle
 931 !      if we have information on this qpt
 932 !      iqpt is equivalent to ieqqpt: get it from file or memory
 933        gkk_qpt_tmp(:,:,:) = gamma_qpt(:,:,:,ieqqpt)
 934 
 935        neqqpt=neqqpt+1
 936 
 937 !
 938 !      MJV note 02/2010:
 939 !      the correspondence of symrel and symrec in the different cases, symmetrizing there
 940 !      and back, has been fixed in the cases with and without scalprod (ie cartesian
 941 !      and reduced real space coordinates) with respect to a calculation with no symmetries
 942 !      I believe everything is settled, but still do not know why the 2 versions of the ss
 943 !      matrices here use different rel/rec, instead of just being multiplied by the rprim gprim...
 944 !
 945        if (ep_scalprod==1) then
 946          do ii=1,3
 947            do jj=1,3
 948              ss(ii,jj)=zero
 949              do kk=1,3
 950                do ll=1,3
 951                  ss(ii,jj)=ss(ii,jj)+Cryst%rprimd(ii,kk)*Cryst%symrel(kk,ll,isym)*Cryst%gprimd(ll,jj)
 952                end do
 953              end do
 954            end do
 955          end do
 956        else
 957          do ii=1,3
 958            do jj=1,3
 959              ss(ii,jj) = Cryst%symrec(ii,jj,isym)
 960            end do
 961          end do
 962        end if
 963 
 964        ss_allatoms(:,:,:) = zero
 965        do iatom=1,natom
 966          ancestor_iatom = Cryst%indsym(4,isym,iatom)
 967          ss_allatoms(1, (ancestor_iatom-1)*3+1:(ancestor_iatom-1)*3+3, (iatom-1)*3+1:(iatom-1)*3+3) = ss(1:3,1:3)
 968        end do
 969 
 970 
 971 !      NOTE   ssinv(ii,jj)=ssinv(ii,jj)+Cryst%gprimd(ii,kk)*rprimd(jj,ll)*Cryst%symrec(ll,kk,isym)
 972 
 973        do isppol=1,nsppol
 974 !        multiply by the ss matrices
 975          tmp_mat2(:,:,:) = zero
 976          tmp_mat(:,:,:) = reshape(gkk_qpt_tmp(:,:,isppol),(/2,nbranch,nbranch/))
 977 
 978          call ZGEMM ('N','N',nbranch,nbranch,nbranch,&
 979 &         c_one,ss_allatoms,nbranch,tmp_mat,nbranch,c_zero,tmp_mat2,nbranch)
 980 
 981          call ZGEMM ('N','T',nbranch,nbranch,nbranch,&
 982 &         c_one,tmp_mat2,nbranch,ss_allatoms,nbranch,c_zero,tmp_mat,nbranch)
 983 
 984 !        add to gkk_qpt_new
 985          do ibranch =1,nbranch
 986            do jbranch =1,nbranch
 987              gkk_qpt_new(:,(jbranch-1)*nbranch+ibranch,isppol) = &
 988 &             gkk_qpt_new(:,(jbranch-1)*nbranch+ibranch,isppol) + tmp_mat(:,jbranch,ibranch)
 989            end do
 990          end do
 991        end do ! isppol
 992 !
 993      end do ! isym
 994    end do ! itim
 995 !
 996    ABI_CHECK(neqqpt>0,'no q-points found equivalent to iqpt ')
 997 !  Divide by number of equivalent qpts found.
 998    gkk_qpt_new(:,:,:) = gkk_qpt_new(:,:,:)/neqqpt
 999 
1000 !  copy the symmetrized version into all the equivalent qpoints, appropriately transformed
1001    do itim=1,2
1002      do isym=1,nsym
1003 !      ieqqpt is sent onto iqpt by itim/isym
1004        ieqqpt = qpttoqpt(itim,isym,iqpt)
1005 
1006        if (symmetrized_qpt(ieqqpt) /= -1) cycle
1007        gkk_qpt_tmp(:,:,:) = zero
1008 
1009 !      use symrec matrices to get inverse transform from isym^{-1}
1010        if (ep_scalprod==1) then
1011          do ii=1,3
1012            do jj=1,3
1013              ss(ii,jj)=zero
1014              do kk=1,3
1015                do ll=1,3
1016 !                Use inverse of symop matrix here to get back to ieqqpt (inv+transpose is in symrec and in gprimd)
1017                  ss(ii,jj)=ss(ii,jj)+Cryst%rprimd(ii,kk)*Cryst%symrec(ll,kk,isym)*Cryst%gprimd(ll,jj)
1018                end do
1019              end do
1020            end do
1021          end do
1022        else
1023          do ii=1,3
1024            do jj=1,3
1025              ss(ii,jj) = Cryst%symrel(jj,ii,isym)
1026            end do
1027          end do
1028        end if
1029 
1030        ss_allatoms(:,:,:) = zero
1031        do iatom=1,natom
1032          ancestor_iatom = Cryst%indsym(4,isym,iatom)
1033          ss_allatoms(1, (ancestor_iatom-1)*3+1:(ancestor_iatom-1)*3+3, (iatom-1)*3+1:(iatom-1)*3+3) = ss(1:3,1:3)
1034        end do
1035 
1036 !      ! Use inverse of symop matrix here to get back to ieqqpt
1037 !      ssinv(ii,jj)=ssinv(ii,jj)+gprimd(ii,kk)*rprimd(jj,ll)*Cryst%symrel(kk,ll,isym)
1038 
1039        do isppol=1,nsppol
1040 !        multiply by the ss^{-1} matrices
1041          tmp_mat2(:,:,:) = zero
1042          tmp_mat(:,:,:) = reshape(gkk_qpt_new(:,:,isppol),(/2,nbranch,nbranch/))
1043 
1044          call ZGEMM ('N','N',nbranch,nbranch,nbranch,&
1045 &         c_one,ss_allatoms,nbranch,tmp_mat,nbranch,c_zero,tmp_mat2,nbranch)
1046 
1047          call ZGEMM ('N','T',nbranch,nbranch,nbranch,&
1048 &         c_one,tmp_mat2,nbranch,ss_allatoms,nbranch,c_zero,tmp_mat,nbranch)
1049 
1050 !        FIXME: the following could just be a reshape
1051          do ibranch =1,nbranch
1052            do jbranch =1,nbranch
1053              gkk_qpt_tmp(:,(jbranch-1)*nbranch+ibranch,isppol) =&
1054 &             tmp_mat(:,jbranch,ibranch)
1055            end do
1056          end do
1057          if (gkk_flag (1,1,isppol,ieqqpt) == -1) gkk_flag (:,:,isppol,ieqqpt) = 0
1058        end do ! end isppol do
1059 
1060 !      save symmetrized matrices for qpt ieqqpt
1061        gamma_qpt(:,:,:,ieqqpt) = gkk_qpt_tmp(:,:,:)
1062 
1063        symmetrized_qpt(ieqqpt) = 1
1064 
1065      end do !isym
1066    end do !itim
1067  end do !iqpt
1068 
1069  ABI_DEALLOCATE(gkk_qpt_new)
1070  ABI_DEALLOCATE(gkk_qpt_tmp)
1071 
1072 end subroutine complete_gamma

ABINIT/complete_gamma_tr [ Functions ]

[ Top ] [ Functions ]

NAME

 complete_gamma_tr

FUNCTION

 Use the set of special q points calculated by the Monkhorst & Pack Technique.
 Check if all the informations for the q points are present in
 the input gamma transport matrices.
 Generate the gamma transport matrices (already summed over the FS) of the set of q points which
 samples homogeneously the entire Brillouin zone.

INPUTS

 crystal<crystal_t>=data type gathering info on the crystalline structure.
 ep_scalprod= flag for scalar product of gkk with phonon displacement vectors
 nbranch=number of phonon branches = 3*natom
 nqptirred=nqpt irred BZ
 nqpt_full=nqpt full BZ
 nsppol=number of spins
 qirredtofull= mapping irred to full qpoints
 qpttoqpt = qpoint index mapping under symops

OUTPUT

 gamma_qpt_tr = in/out: set of gamma matrix elements completed and symmetrized
    gamma_qpt_tr(2,9,nbranch*nbranch,nsppol,nqpt_full)

PARENTS

      elphon

CHILDREN

      dgemm

SOURCE

1109 subroutine complete_gamma_tr(crystal,ep_scalprod,nbranch,nqptirred,nqpt_full,nsppol,gamma_qpt_tr,qirredtofull,qpttoqpt)
1110 
1111  use m_linalg_interfaces
1112 
1113 !This section has been created automatically by the script Abilint (TD).
1114 !Do not modify the following lines by hand.
1115 #undef ABI_FUNC
1116 #define ABI_FUNC 'complete_gamma_tr'
1117 !End of the abilint section
1118 
1119  implicit none
1120 
1121 !Arguments ------------------------------------
1122 !scalars
1123  integer, intent(in) :: nbranch,nqptirred,nqpt_full,nsppol, ep_scalprod
1124  type(crystal_t),intent(in) :: crystal
1125 !arrays
1126  integer,intent(in) :: qpttoqpt(2,crystal%nsym,nqpt_full)
1127  integer,intent(in) :: qirredtofull(nqptirred)
1128  real(dp), intent(inout) :: gamma_qpt_tr(2,9,nbranch*nbranch,nsppol,nqpt_full)
1129 
1130 !Local variables-------------------------------
1131 !scalars
1132  integer :: ieqqpt,ii,iqpt,isppol,isym
1133  integer :: itim,jj,kk,ll,neqqpt
1134  integer :: iatom,ancestor_iatom
1135  integer :: iqpt_fullbz,imode, itensor
1136  real(dp),parameter :: tol=2.d-8
1137 !arrays
1138  integer :: symrel(3,3,crystal%nsym),symrec(3,3,crystal%nsym)
1139  integer :: symmetrized_qpt(nqpt_full)
1140  integer :: gkk_flag(nbranch,nbranch,nsppol,nqpt_full)
1141  real(dp) :: gprimd(3,3),rprimd(3,3)
1142  real(dp) :: ss(3,3), sscart(3,3)
1143  real(dp) :: tmp_mat(2,nbranch,nbranch)
1144  real(dp) :: tmp_mat2(2,nbranch,nbranch)
1145  real(dp) :: tmp_tensor(2,3,3)
1146  real(dp) :: tmp_tensor2(2,3,3)
1147  real(dp) :: ss_allatoms(nbranch,nbranch)
1148  real(dp) :: c_one(2), c_zero(2)
1149  real(dp),allocatable :: gkk_qpt_new(:,:,:,:),gkk_qpt_tmp(:,:,:,:)
1150 
1151 ! *********************************************************************
1152 
1153  c_one = (/one,zero/)
1154  c_zero = (/zero,zero/)
1155 
1156  gprimd = crystal%gprimd
1157  rprimd = crystal%rprimd
1158 
1159  symrec =  crystal%symrec
1160  symrel =  crystal%symrel
1161 
1162 !Generation of the gkk matrices relative to the q points
1163 !of the set which samples the entire Brillouin zone
1164 
1165 !set up flags for gamma_qpt matrices we have
1166  gkk_flag = -1
1167  do iqpt=1,nqptirred
1168    iqpt_fullbz = qirredtofull(iqpt)
1169    gkk_flag(:,:,:,iqpt_fullbz) = 1
1170  end do
1171 
1172  symmetrized_qpt(:) = -1
1173 ! isppol=1
1174 
1175  ABI_ALLOCATE(gkk_qpt_new,(2,9,nbranch*nbranch, nsppol))
1176  ABI_ALLOCATE(gkk_qpt_tmp,(2,9,nbranch*nbranch, nsppol))
1177 
1178  do iqpt=1,nqpt_full
1179 
1180 !  Already symmetrized?
1181    if (symmetrized_qpt(iqpt) == 1) cycle
1182 
1183    gkk_qpt_new(:,:,:,:) = zero
1184 
1185 !  loop over qpoints equivalent to iqpt
1186    neqqpt=0
1187 !  do not use time reversal symmetry to complete the qpoints:
1188 !  do not know what happens to the gamma matrices
1189 
1190    do itim=1,2
1191      do isym=1,crystal%nsym
1192 !      ieqqpt is sent onto iqpt by itim/isym
1193        ieqqpt = qpttoqpt(itim,isym,iqpt)
1194 
1195        if (gkk_flag(1,1,1,ieqqpt) == -1) cycle
1196 !      if we have information on this qpt
1197 !      iqpt is equivalent to ieqqpt: get it from file or memory
1198        gkk_qpt_tmp(:,:,:,:) = gamma_qpt_tr(:,:,:,:,ieqqpt)
1199 
1200        neqqpt=neqqpt+1
1201 
1202 !
1203 !      MJV note 02/2010:
1204 !      the correspondence of symrel and symrec in the different cases, symmetrizing there
1205 !      and back, has been fixed in the cases with and without scalprod (ie cartesian
1206 !      and reduced real space coordinates) with respect to a calculation with no symmetries
1207 !      I believe everything is settled, but still do not know why the 2 versions of the ss
1208 !      matrices here use different rel/rec, instead of just being multiplied by the rprim gprim...
1209 !
1210        do ii=1,3
1211          do jj=1,3
1212            sscart(ii,jj)=0.0_dp
1213            do kk=1,3
1214              do ll=1,3
1215                sscart(ii,jj)=sscart(ii,jj)+rprimd(ii,kk)*symrel(kk,ll,isym)*gprimd(ll,jj)
1216 !              sscart(ii,jj)=sscart(ii,jj)+rprimd(ii,kk)*symrel(kk,ll,isym)*gprimd(ll,jj)
1217              end do
1218            end do
1219          end do
1220        end do
1221        if (ep_scalprod==1) then
1222          do ii=1,3
1223            do jj=1,3
1224              ss(ii,jj)=0.0_dp
1225              do kk=1,3
1226                do ll=1,3
1227                  ss(ii,jj)=ss(ii,jj)+rprimd(ii,kk)*symrel(kk,ll,isym)*gprimd(ll,jj)
1228                end do
1229              end do
1230            end do
1231          end do
1232        else
1233          do ii=1,3
1234            do jj=1,3
1235              ss(ii,jj) = symrec(ii,jj,isym)
1236            end do
1237          end do
1238        end if
1239 
1240        ss_allatoms(:,:) = zero
1241        do iatom=1,crystal%natom
1242          ancestor_iatom = crystal%indsym(4,isym,iatom)
1243          ss_allatoms((ancestor_iatom-1)*3+1:(ancestor_iatom-1)*3+3,&
1244 &         (iatom-1)*3+1:         (iatom-1)*3+3) = ss(1:3,1:3)
1245        end do
1246 
1247 
1248 !      NOTE   ssinv(ii,jj)=ssinv(ii,jj)+gprimd(ii,kk)*rprimd(jj,ll)*symrec(ll,kk,isym)
1249 
1250        do isppol=1,nsppol
1251 
1252 !        for each tensor component, rotate the cartesian directions of phonon modes
1253          do itensor = 1, 9
1254 !          multiply by the ss matrices
1255            tmp_mat2(:,:,:) = zero
1256            tmp_mat(:,:,:) = reshape(gkk_qpt_tmp(:,itensor,:,isppol),&
1257 &           (/2,nbranch,nbranch/))
1258            call DGEMM ('N','N',nbranch,nbranch,nbranch,&
1259 &           one,ss_allatoms,nbranch,tmp_mat(1,:,:),nbranch,zero,&
1260 &           tmp_mat2(1,:,:),nbranch)
1261            call DGEMM ('N','N',nbranch,nbranch,nbranch,&
1262 &           one,ss_allatoms,nbranch,tmp_mat(2,:,:),nbranch,zero,&
1263 &           tmp_mat2(2,:,:),nbranch)
1264 
1265            call DGEMM ('N','T',nbranch,nbranch,nbranch,&
1266 &           one,tmp_mat2(1,:,:),nbranch,ss_allatoms,nbranch,zero,&
1267 &           tmp_mat(1,:,:),nbranch)
1268            call DGEMM ('N','T',nbranch,nbranch,nbranch,&
1269 &           one,tmp_mat2(2,:,:),nbranch,ss_allatoms,nbranch,zero,&
1270 &           tmp_mat(2,:,:),nbranch)
1271 
1272            gkk_qpt_tmp(:,itensor,:,isppol) = reshape (tmp_mat, (/2,nbranch*nbranch/))
1273          end do ! itensor
1274 
1275 !        for each cartesian direction/phonon mode, rotate the tensor components
1276          do imode = 1, nbranch*nbranch
1277            tmp_tensor2(:,:,:) = zero
1278            tmp_tensor(:,:,:) = reshape(gkk_qpt_tmp(:,:,imode,isppol),&
1279 &           (/2,3,3/))
1280            call DGEMM ('N','N',3,3,3,&
1281 &           one,sscart,3,tmp_tensor(1,:,:),3,zero,&
1282 &           tmp_tensor2(1,:,:),3)
1283            call DGEMM ('N','T',3,3,3,&
1284 &           one,tmp_tensor2(1,:,:),3,sscart,3,zero,&
1285 &           tmp_tensor(1,:,:),3)
1286 
1287            call DGEMM ('N','N',3,3,3,&
1288 &           one,sscart,3,tmp_tensor(2,:,:),3,zero,&
1289 &           tmp_tensor2(2,:,:),3)
1290            call DGEMM ('N','T',3,3,3,&
1291 &           one,tmp_tensor2(2,:,:),3,sscart,3,zero,&
1292 &           tmp_tensor(2,:,:),3)
1293 
1294            gkk_qpt_tmp(:,:,imode,isppol) = reshape (tmp_tensor, (/2,9/)) ! modified by BX
1295          end do ! imode
1296 
1297 !        add to gkk_qpt_new
1298          gkk_qpt_new(:,:,:,isppol) = gkk_qpt_new(:,:,:,isppol) + gkk_qpt_tmp(:,:,:,isppol)
1299 
1300        end do ! end isppol do
1301 
1302      end do ! end isym do
1303    end do ! end itim do
1304 
1305    ABI_CHECK(neqqpt>0,'no q-points found equivalent to iqpt ')
1306 
1307 !  divide by number of equivalent qpts found
1308    gkk_qpt_new = gkk_qpt_new/neqqpt
1309 
1310 
1311 !  copy the symmetrized version into all the equivalent qpoints, appropriately transformed
1312    do itim=1,2
1313      do isym=1,crystal%nsym
1314 !      ieqqpt is sent onto iqpt by itim/isym
1315        ieqqpt = qpttoqpt(itim,isym,iqpt)
1316 
1317        if (symmetrized_qpt(ieqqpt) /= -1) cycle
1318        gkk_qpt_tmp = zero
1319 
1320 !      use symrec matrices to get inverse transform from isym^{-1}
1321        do ii=1,3
1322          do jj=1,3
1323            sscart(ii,jj)=0.0_dp
1324            do kk=1,3
1325              do ll=1,3
1326 !              Use inverse of symop matrix here to get back to ieqqpt (inv+transpose is in symrec and in gprimd)
1327 !              sscart(ii,jj)=sscart(ii,jj)+rprimd(ii,kk)*symrec(ll,kk,isym)*gprimd(ll,jj)
1328                sscart(ii,jj)=sscart(ii,jj)+rprimd(ii,kk)*symrec(ll,kk,isym)*gprimd(ll,jj)
1329              end do
1330            end do
1331          end do
1332        end do
1333        if (ep_scalprod==1) then
1334          do ii=1,3
1335            do jj=1,3
1336              ss(ii,jj)=0.0_dp
1337              do kk=1,3
1338                do ll=1,3
1339 !                Use inverse of symop matrix here to get back to ieqqpt (inv+transpose is in symrec and in gprimd)
1340                  ss(ii,jj)=ss(ii,jj)+rprimd(ii,kk)*symrec(ll,kk,isym)*gprimd(ll,jj)
1341                end do
1342              end do
1343            end do
1344          end do
1345        else
1346          do ii=1,3
1347            do jj=1,3
1348              ss(ii,jj) = symrel(jj,ii,isym)
1349            end do
1350          end do
1351        end if
1352 
1353        ss_allatoms(:,:) = zero
1354        do iatom=1,crystal%natom
1355          ancestor_iatom = crystal%indsym(4,isym,iatom)
1356          ss_allatoms((ancestor_iatom-1)*3+1:(ancestor_iatom-1)*3+3,&
1357 &         (iatom-1)*3+1:          (iatom-1)*3+3) = ss(1:3,1:3)
1358        end do
1359 
1360 !      ! Use inverse of symop matrix here to get back to ieqqpt
1361 !      ssinv(ii,jj)=ssinv(ii,jj)+gprimd(ii,kk)*rprimd(jj,ll)*symrel(kk,ll,isym)
1362 
1363        do isppol=1,nsppol
1364          do itensor = 1, 9
1365 !          multiply by the ss^{-1} matrices
1366            tmp_mat2(:,:,:) = zero
1367            tmp_mat(:,:,:) = reshape(gkk_qpt_new(:,itensor,:,isppol),&
1368 &           (/2,nbranch,nbranch/))
1369 
1370 
1371            call DGEMM ('N','N',nbranch,nbranch,nbranch,&
1372 &           one,ss_allatoms,nbranch,tmp_mat(1,:,:),nbranch,zero,&
1373 &           tmp_mat2(1,:,:),nbranch)
1374            call DGEMM ('N','N',nbranch,nbranch,nbranch,&
1375 &           one,ss_allatoms,nbranch,tmp_mat(2,:,:),nbranch,zero,&
1376 &           tmp_mat2(2,:,:),nbranch)
1377 
1378            call DGEMM ('N','T',nbranch,nbranch,nbranch,&
1379 &           one,tmp_mat2(1,:,:),nbranch,ss_allatoms,nbranch,zero,&
1380 &           tmp_mat(1,:,:),nbranch)
1381            call DGEMM ('N','T',nbranch,nbranch,nbranch,&
1382 &           one,tmp_mat2(2,:,:),nbranch,ss_allatoms,nbranch,zero,&
1383 &           tmp_mat(2,:,:),nbranch)
1384 
1385 
1386            gkk_qpt_tmp(:,itensor,:,isppol) = reshape (tmp_mat, (/2,nbranch*nbranch/))
1387          end do ! itensor
1388 
1389 !        for each cartesian direction/phonon mode, rotate the tensor components
1390          do imode = 1, nbranch*nbranch
1391            tmp_tensor2(:,:,:) = zero
1392            tmp_tensor(:,:,:) = reshape(gkk_qpt_tmp(:,:,imode,isppol),&
1393 &           (/2,3,3/))
1394            call DGEMM ('N','N',3,3,3,&
1395 &           one,sscart,3,tmp_tensor(1,:,:),3,zero,&
1396 &           tmp_tensor2(1,:,:),3)
1397            call DGEMM ('N','T',3,3,3,&
1398 &           one,tmp_tensor2(1,:,:),3,sscart,3,zero,&
1399 &           tmp_tensor(1,:,:),3)
1400 
1401            call DGEMM ('N','N',3,3,3,&
1402 &           one,sscart,3,tmp_tensor(2,:,:),3,zero,&
1403 &           tmp_tensor2(2,:,:),3)
1404            call DGEMM ('N','T',3,3,3,&
1405 &           one,tmp_tensor2(2,:,:),3,sscart,3,zero,&
1406 &           tmp_tensor(2,:,:),3)
1407 
1408 !          gkk_qpt_new(:,:,imode,isppol) = reshape (tmp_tensor, (/2,9/)) ! Modified by BX
1409            gkk_qpt_tmp(:,:,imode,isppol) = reshape (tmp_tensor, (/2,9/)) ! Modified by BX
1410          end do ! imode
1411 
1412          if (gkk_flag (1,1,isppol,ieqqpt) == -1) then
1413            gkk_flag (:,:,isppol,ieqqpt) = 0
1414          end if
1415 
1416        end do ! end isppol do
1417 
1418 
1419 !      save symmetrized matrices for qpt ieqqpt
1420        gamma_qpt_tr(:,:,:,:,ieqqpt) = gkk_qpt_tmp(:,:,:,:)
1421 
1422        symmetrized_qpt(ieqqpt) = 1
1423 
1424      end do ! end isym do
1425    end do ! end itim do
1426 
1427  end do
1428 !end iqpt do
1429 
1430  ABI_DEALLOCATE(gkk_qpt_new)
1431  ABI_DEALLOCATE(gkk_qpt_tmp)
1432 
1433 end subroutine complete_gamma_tr

ABINIT/defs_elphon [ Modules ]

[ Top ] [ Modules ]

NAME

 defs_elphon

FUNCTION

 This module contains the datastructures for elphon
  the different (huge) matrices will either be allocated and
  used, or be written to disk. All combinations should be feasible.

COPYRIGHT

 Copyright (C) 2004-2018 ABINIT group (MVer, MG)
 This file is distributed under the terms of the
 GNU General Public Licence, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

NOTES

  Contains the following datastructures:
   1) elph_type contains data and dimensions for the kpoints near the
      fermi surface and the $g_{k k+q}$ matrix elements

PARENTS

CHILDREN

SOURCE

30 #if defined HAVE_CONFIG_H
31 #include "config.h"
32 #endif
33 
34 #include "abi_common.h"
35 
36 module defs_elphon
37 
38  use defs_basis
39  use m_abicore
40  use m_errors
41  use m_xmpi
42 
43  use m_kptrank,  only : kptrank_type, destroy_kptrank, copy_kptrank
44  use m_crystal,    only : crystal_t
45 
46  implicit none
47 
48  private

defs_elphon/elph_ds_clean [ Functions ]

[ Top ] [ defs_elphon ] [ Functions ]

NAME

   elph_ds_clean

FUNCTION

   deallocate remaining arrays in the elph_ds datastructure

INPUTS

  elph_ds = elphon datastructure

PARENTS

      elphon

CHILDREN

NOTES

SOURCE

323 subroutine elph_ds_clean(elph_ds)
324 
325 
326 !This section has been created automatically by the script Abilint (TD).
327 !Do not modify the following lines by hand.
328 #undef ABI_FUNC
329 #define ABI_FUNC 'elph_ds_clean'
330 !End of the abilint section
331 
332  implicit none
333 
334 !Arguments ------------------------------------
335 !scalars
336  type(elph_type), intent(inout) :: elph_ds
337 
338 ! *************************************************************************
339 
340  !@elph_type
341  if (allocated(elph_ds%qirredtofull))  then
342    ABI_DEALLOCATE(elph_ds%qirredtofull)
343  end if
344  if (allocated(elph_ds%wtq))  then
345    ABI_DEALLOCATE(elph_ds%wtq)
346  end if
347  if (allocated(elph_ds%n0))  then
348    ABI_DEALLOCATE(elph_ds%n0)
349  end if
350  if (allocated(elph_ds%qpt_full))  then
351    ABI_DEALLOCATE(elph_ds%qpt_full)
352  end if
353  if (allocated(elph_ds%gkk_intweight))  then
354    ABI_DEALLOCATE(elph_ds%gkk_intweight)
355  end if
356  if (allocated(elph_ds%gkk_qpt))  then
357    ABI_DEALLOCATE(elph_ds%gkk_qpt)
358  end if
359  if (allocated(elph_ds%gkk_rpt))  then
360    ABI_DEALLOCATE(elph_ds%gkk_rpt)
361  end if
362  if (allocated(elph_ds%gkk2))  then
363    ABI_DEALLOCATE(elph_ds%gkk2)
364  end if
365  if (allocated(elph_ds%gamma_qpt))  then
366    ABI_DEALLOCATE(elph_ds%gamma_qpt)
367  end if
368  if (allocated(elph_ds%gamma_rpt))  then
369    ABI_DEALLOCATE(elph_ds%gamma_rpt)
370  end if
371  if (allocated(elph_ds%phfrq))  then
372    ABI_DEALLOCATE(elph_ds%phfrq)
373  end if
374  if (allocated(elph_ds%a2f))  then
375    ABI_DEALLOCATE(elph_ds%a2f)
376  end if
377  if (allocated(elph_ds%qgrid_data))  then
378    ABI_DEALLOCATE(elph_ds%qgrid_data)
379  end if
380 
381  call elph_k_destroy (elph_ds%k_phon)
382  call elph_k_destroy (elph_ds%k_fine)
383 
384  call destroy_kptrank (elph_ds%k_fine%kptrank_t)
385 
386 end subroutine elph_ds_clean

defs_elphon/elph_k_copy [ Functions ]

[ Top ] [ defs_elphon ] [ Functions ]

NAME

   elph_k_copy

FUNCTION

   allocate and copy arrays in the elph_k datastructure

INPUTS

  elph_k = elphon k-points datastructure

PARENTS

CHILDREN

NOTES

SOURCE

532 subroutine elph_k_copy(elph_k_in, elph_k_out)
533 
534 
535 !This section has been created automatically by the script Abilint (TD).
536 !Do not modify the following lines by hand.
537 #undef ABI_FUNC
538 #define ABI_FUNC 'elph_k_copy'
539 !End of the abilint section
540 
541  implicit none
542 
543 !Arguments ------------------------------------
544 !scalars
545  type(elph_kgrid_type), intent(in) :: elph_k_in
546  type(elph_kgrid_type), intent(out) :: elph_k_out
547 
548 ! *************************************************************************
549 
550  !@elph_kgrid_type
551  elph_k_out%nband = elph_k_in%nband
552  elph_k_out%nsppol = elph_k_in%nsppol
553  elph_k_out%nsym = elph_k_in%nsym
554 
555  elph_k_out%nkpt = elph_k_in%nkpt
556  elph_k_out%nkptirr = elph_k_in%nkptirr
557 
558  elph_k_out%my_nkpt = elph_k_in%my_nkpt
559 
560  ABI_ALLOCATE(elph_k_out%my_kpt,(elph_k_out%nkpt))
561  elph_k_out%my_kpt = elph_k_in%my_kpt
562 
563  ABI_ALLOCATE(elph_k_out%my_ikpt,(elph_k_out%my_nkpt))
564  elph_k_out%my_ikpt = elph_k_in%my_ikpt
565 
566  ABI_ALLOCATE(elph_k_out%kptirr,(3,elph_k_out%nkptirr))
567  elph_k_out%kptirr = elph_k_in%kptirr
568  ABI_ALLOCATE(elph_k_out%wtkirr,(elph_k_out%nkptirr))
569  elph_k_out%wtkirr = elph_k_in%wtkirr
570 
571  ABI_ALLOCATE(elph_k_out%wtk,(elph_k_out%nband,elph_k_out%nkpt,elph_k_out%nsppol))
572  elph_k_out%wtk = elph_k_in%wtk
573  ABI_ALLOCATE(elph_k_out%kpt,(3,elph_k_out%nkpt))
574  elph_k_out%kpt = elph_k_in%kpt
575 
576  call copy_kptrank(elph_k_in%kptrank_t, elph_k_out%kptrank_t)
577 
578  ABI_ALLOCATE(elph_k_out%irr2full,(elph_k_out%nkptirr))
579  elph_k_out%irr2full = elph_k_in%irr2full
580  ABI_ALLOCATE(elph_k_out%full2irr,(3,elph_k_out%nkpt))
581  elph_k_out%full2irr = elph_k_in%full2irr
582  ABI_ALLOCATE(elph_k_out%full2full,(2,elph_k_out%nsym,elph_k_out%nkpt))
583  elph_k_out%full2full = elph_k_in%full2full
584 
585  ABI_ALLOCATE(elph_k_out%irredtoGS,(elph_k_out%nkptirr))
586  elph_k_out%irredtoGS = elph_k_in%irredtoGS
587 
588 
589 end subroutine elph_k_copy

defs_elphon/elph_k_destroy [ Functions ]

[ Top ] [ defs_elphon ] [ Functions ]

NAME

   elph_k_destroy

FUNCTION

   deallocate arrays in the elph_k datastructure

INPUTS

  elph_k = elphon k-points datastructure

PARENTS

      defs_elphon

CHILDREN

NOTES

SOURCE

614 subroutine elph_k_destroy(elph_k)
615 
616 
617 !This section has been created automatically by the script Abilint (TD).
618 !Do not modify the following lines by hand.
619 #undef ABI_FUNC
620 #define ABI_FUNC 'elph_k_destroy'
621 !End of the abilint section
622 
623  implicit none
624 
625 !Arguments ------------------------------------
626 !scalars
627  type(elph_kgrid_type), intent(inout) :: elph_k
628 
629 ! *************************************************************************
630 
631  !@elph_kgrid_type
632  if(allocated(elph_k%irr2full)) then
633    ABI_DEALLOCATE (elph_k%irr2full)
634  end if
635  if(allocated(elph_k%full2irr)) then
636    ABI_DEALLOCATE (elph_k%full2irr)
637  end if
638  if(allocated(elph_k%full2full)) then
639    ABI_DEALLOCATE (elph_k%full2full)
640  end if
641  if(allocated(elph_k%irredtoGS)) then
642    ABI_DEALLOCATE (elph_k%irredtoGS)
643  end if
644  if(allocated(elph_k%new_irredtoGS)) then
645    ABI_DEALLOCATE (elph_k%new_irredtoGS)
646  end if
647  if(allocated(elph_k%kpt)) then
648    ABI_DEALLOCATE (elph_k%kpt)
649  end if
650  if(allocated(elph_k%kptirr)) then
651    ABI_DEALLOCATE (elph_k%kptirr)
652  end if
653  if(allocated(elph_k%new_kptirr)) then
654    ABI_DEALLOCATE (elph_k%new_kptirr)
655  end if
656  if(allocated(elph_k%my_kpt)) then
657    ABI_DEALLOCATE (elph_k%my_kpt)
658  end if
659  if(allocated(elph_k%my_ikpt)) then
660    ABI_DEALLOCATE (elph_k%my_ikpt)
661  end if
662  if(allocated(elph_k%wtk)) then
663    ABI_DEALLOCATE (elph_k%wtk)
664  end if
665  if(allocated(elph_k%wtq)) then
666    ABI_DEALLOCATE (elph_k%wtq)
667  end if
668  if(allocated(elph_k%wtkirr)) then
669    ABI_DEALLOCATE (elph_k%wtkirr)
670  end if
671  if(allocated(elph_k%new_wtkirr)) then
672    ABI_DEALLOCATE (elph_k%new_wtkirr)
673  end if
674  if(allocated(elph_k%velocwtk)) then
675    ABI_DEALLOCATE (elph_k%velocwtk)
676  end if
677  if(allocated(elph_k%vvelocwtk)) then
678    ABI_DEALLOCATE (elph_k%vvelocwtk)
679  end if
680 
681  call destroy_kptrank (elph_k%kptrank_t)
682 
683 end subroutine elph_k_destroy

defs_elphon/elph_k_procs [ Functions ]

[ Top ] [ defs_elphon ] [ Functions ]

NAME

   elph_k_procs

FUNCTION

   allocate kpt to processors, in the elph_k datastructure

INPUTS

  nproc = number of k-parallel processors
  elph_k = elphon k-points datastructure

PARENTS

      elphon

CHILDREN

NOTES

SOURCE

709 subroutine elph_k_procs(nproc, elph_k)
710 
711 
712 !This section has been created automatically by the script Abilint (TD).
713 !Do not modify the following lines by hand.
714 #undef ABI_FUNC
715 #define ABI_FUNC 'elph_k_procs'
716 !End of the abilint section
717 
718  implicit none
719 
720 !Arguments ------------------------------------
721 !scalars
722  integer, intent(in) :: nproc
723  type(elph_kgrid_type), intent(inout) :: elph_k
724 
725  integer :: ikpt, me, ik_this_proc
726 ! *************************************************************************
727 
728  if (allocated(elph_k%my_kpt)) then
729    ABI_DEALLOCATE (elph_k%my_kpt)
730  end if
731  ABI_ALLOCATE (elph_k%my_kpt, (elph_k%nkpt))
732 
733  elph_k%my_kpt = 0
734  elph_k%my_nkpt = 0
735  me = xmpi_comm_rank(xmpi_world)
736  do ikpt = 1, elph_k%nkpt
737    elph_k%my_kpt(ikpt) = MOD(ikpt-1, nproc)
738    if (elph_k%my_kpt(ikpt) == me) elph_k%my_nkpt = elph_k%my_nkpt + 1
739  end do
740 
741 ! create inverse mapping from ik_this_proc to ikpt
742  if (allocated(elph_k%my_ikpt)) then
743    ABI_DEALLOCATE (elph_k%my_ikpt)
744  end if
745  ABI_ALLOCATE (elph_k%my_ikpt, (elph_k%my_nkpt))
746  elph_k%my_ikpt = 0
747 
748  ik_this_proc = 0
749  do ikpt = 1, elph_k%nkpt
750    if (elph_k%my_kpt(ikpt) == me) then
751      ik_this_proc = ik_this_proc + 1
752      elph_k%my_ikpt(ik_this_proc) = ikpt
753    end if
754  end do
755  ABI_CHECK(ik_this_proc == elph_k%my_nkpt, 'found inconsistent k distribution in processors')
756 
757  write (std_out,*) 'elph_k_procs : nkpt, distrib = ', elph_k%my_nkpt
758  write (std_out,*) elph_k%my_kpt
759  write (std_out,*) elph_k%my_ikpt
760 
761 end subroutine elph_k_procs

defs_elphon/elph_kgrid_type [ Types ]

[ Top ] [ defs_elphon ] [ Types ]

NAME

 elph_kgrid_type

FUNCTION

 elph_kgrid_type contains k-point grid data and dimensions
  this is a sub object of elph_type

SOURCE

 66   type,public :: elph_kgrid_type
 67 
 68    integer :: nband                           ! number of bands for weights
 69    integer :: nsppol                          ! number of spin pol for weights
 70    integer :: nsym                            ! number of symmetry operations
 71    integer :: nkpt                            ! number of k-points in full grid
 72    integer :: nkptirr                         ! number of k-points in irreducible grid
 73    integer :: new_nkptirr                     ! number of k-points in irreducible grid
 74    integer :: my_nkpt                         ! number of k-points on present processor
 75 
 76    type(kptrank_type) :: kptrank_t            ! ranking of all kpoints on phonon calculation grid, and inverse rank
 77 
 78    integer, allocatable :: irr2full(:)            ! correspondence of irred kpoints to a full one
 79    integer, allocatable :: full2irr(:,:)          ! correspondence of full k to one irred kpoints through sym and timrev
 80    integer, allocatable :: full2full(:,:,:)       ! symmetry mapping of kpoints
 81    integer, allocatable :: my_kpt(:)              ! flag for k-points belonging to present proc (= me index of proc for each k-point)
 82    integer, allocatable :: my_ikpt(:)             ! flag for k-points belonging to present proc (= me index of proc for each k-point)
 83    integer, allocatable :: irredtoGS(:)           ! (nkptirr)
 84    integer, allocatable :: new_irredtoGS(:)       ! (new_nkptirr)
 85 
 86    real(dp),allocatable :: kpt(:,:)               ! coordinates of the full kpoints from phonon calculation
 87    real(dp),allocatable :: kptirr(:,:)            ! irreducible k-points, for preliminary set up
 88    real(dp),allocatable :: new_kptirr(:,:)        ! irreducible k-points, for preliminary set up
 89    real(dp),allocatable :: wtk(:,:,:)             ! integration weights (see also gkk_intweight)
 90    real(dp),allocatable :: wtq(:,:,:)             ! integration weights (see also gkk_intweight)
 91    real(dp),allocatable :: wtkirr(:)              ! weights for irreducible kpoints, to sum over _whole_ BZ (not just Fermi Surface)
 92    real(dp),allocatable :: new_wtkirr(:)          ! weights for irreducible kpoints, to sum over _whole_ BZ (not just Fermi Surface)
 93    real(dp),allocatable :: velocwtk(:,:,:,:)      ! (nFSband,nkpt_fine,3,nsppol), v(k)*wtk
 94    real(dp),allocatable :: vvelocwtk(:,:,:,:,:)   ! (nFSband,nkpt_fine,3,3,nsppol), v(k)*v(k)*wtk
 95 
 96   end type elph_kgrid_type
 97 
 98   public :: elph_k_copy
 99   public :: elph_k_procs
100   public :: elph_k_destroy

defs_elphon/elph_tr_ds_clean [ Functions ]

[ Top ] [ defs_elphon ] [ Functions ]

NAME

   elph_tr_ds_clean

FUNCTION

   deallocate remaining arrays in the elph_tr_ds datastructure

INPUTS

  elph_tr_ds = elphon transport datastructure

PARENTS

      elphon

CHILDREN

NOTES

SOURCE

411 subroutine elph_tr_ds_clean(elph_tr_ds)
412 
413 
414 !This section has been created automatically by the script Abilint (TD).
415 !Do not modify the following lines by hand.
416 #undef ABI_FUNC
417 #define ABI_FUNC 'elph_tr_ds_clean'
418 !End of the abilint section
419 
420  implicit none
421 
422 !Arguments ------------------------------------
423 !scalars
424  type(elph_tr_type), intent(inout) :: elph_tr_ds
425 
426 ! *************************************************************************
427 
428  !@elph_tr_type
429  if (allocated(elph_tr_ds%el_veloc))  then
430    ABI_DEALLOCATE(elph_tr_ds%el_veloc)
431  end if
432  if (allocated(elph_tr_ds%FSelecveloc_sq))  then
433    ABI_DEALLOCATE(elph_tr_ds%FSelecveloc_sq)
434  end if
435  if (allocated(elph_tr_ds%veloc_sq0))  then
436    ABI_DEALLOCATE(elph_tr_ds%veloc_sq0)
437  end if
438  if (allocated(elph_tr_ds%veloc_sq))  then
439    ABI_DEALLOCATE(elph_tr_ds%veloc_sq)
440  end if
441  if (allocated(elph_tr_ds%dos_n0))  then
442    ABI_DEALLOCATE(elph_tr_ds%dos_n0)
443  end if
444  if (allocated(elph_tr_ds%dos_n))  then
445    ABI_DEALLOCATE(elph_tr_ds%dos_n)
446  end if
447  if (allocated(elph_tr_ds%en_all))  then
448    ABI_DEALLOCATE(elph_tr_ds%en_all)
449  end if
450  if (allocated(elph_tr_ds%de_all))  then
451    ABI_DEALLOCATE(elph_tr_ds%de_all)
452  end if
453  if (allocated(elph_tr_ds%gamma_qpt_tr))  then
454    ABI_DEALLOCATE(elph_tr_ds%gamma_qpt_tr)
455  end if
456  if (allocated(elph_tr_ds%gamma_qpt_trin))  then
457    ABI_DEALLOCATE(elph_tr_ds%gamma_qpt_trin)
458  end if
459  if (allocated(elph_tr_ds%gamma_qpt_trout))  then
460    ABI_DEALLOCATE(elph_tr_ds%gamma_qpt_trout)
461  end if
462  if (allocated(elph_tr_ds%gamma_rpt_tr))  then
463    ABI_DEALLOCATE(elph_tr_ds%gamma_rpt_tr)
464  end if
465  if (allocated(elph_tr_ds%gamma_rpt_trin))  then
466    ABI_DEALLOCATE(elph_tr_ds%gamma_rpt_trin)
467  end if
468  if (allocated(elph_tr_ds%gamma_rpt_trout))  then
469    ABI_DEALLOCATE(elph_tr_ds%gamma_rpt_trout)
470  end if
471  if (allocated(elph_tr_ds%a2f_1d_tr))  then
472    ABI_DEALLOCATE(elph_tr_ds%a2f_1d_tr)
473  end if
474  if (allocated(elph_tr_ds%a2f_1d_trin))  then
475    ABI_DEALLOCATE(elph_tr_ds%a2f_1d_trin)
476  end if
477  if (allocated(elph_tr_ds%a2f_1d_trout))  then
478    ABI_DEALLOCATE(elph_tr_ds%a2f_1d_trout)
479  end if
480  if (allocated(elph_tr_ds%tmp_gkk_intweight))  then
481    ABI_DEALLOCATE(elph_tr_ds%tmp_gkk_intweight)
482  end if
483  if (allocated(elph_tr_ds%tmp_gkk_intweight1))  then
484    ABI_DEALLOCATE(elph_tr_ds%tmp_gkk_intweight1)
485  end if
486  if (allocated(elph_tr_ds%tmp_gkk_intweight2))  then
487    ABI_DEALLOCATE(elph_tr_ds%tmp_gkk_intweight2)
488  end if
489  if (allocated(elph_tr_ds%tmp_velocwtk))  then
490    ABI_DEALLOCATE(elph_tr_ds%tmp_velocwtk)
491  end if
492  if (allocated(elph_tr_ds%tmp_velocwtk1))  then
493    ABI_DEALLOCATE(elph_tr_ds%tmp_velocwtk1)
494  end if
495  if (allocated(elph_tr_ds%tmp_velocwtk2))  then
496    ABI_DEALLOCATE(elph_tr_ds%tmp_velocwtk2)
497  end if
498  if (allocated(elph_tr_ds%tmp_vvelocwtk))  then
499    ABI_DEALLOCATE(elph_tr_ds%tmp_vvelocwtk)
500  end if
501  if (allocated(elph_tr_ds%tmp_vvelocwtk1))  then
502    ABI_DEALLOCATE(elph_tr_ds%tmp_vvelocwtk1)
503  end if
504  if (allocated(elph_tr_ds%tmp_vvelocwtk2))  then
505    ABI_DEALLOCATE(elph_tr_ds%tmp_vvelocwtk2)
506  end if
507 
508 end subroutine elph_tr_ds_clean

defs_elphon/elph_tr_type [ Types ]

[ Top ] [ defs_elphon ] [ Types ]

NAME

 elph_tr_type

FUNCTION

 elph_tr_ds contains the necessary data for the transport properties

SOURCE

247   type,public :: elph_tr_type
248 
249      integer :: ifltransport
250      integer :: unitgkq_trin,unitgkq_trout
251      integer :: gkqwrite,gkqexist
252      integer :: onegkksize
253 
254      character(len=fnlen) :: ddkfilename
255 
256      real(dp),allocatable :: dos_n0(:,:)                  ! (nT,nsppol) DOS at the Fermi level (states/Ha/spin) at input temperatures
257      real(dp),allocatable :: dos_n(:,:)                   ! (nE,nsppol) DOS at the selected energies (states/Ha/spin)
258      real(dp),allocatable :: en_all(:,:)                  ! (nE,nsppol) selected energies
259      real(dp),allocatable :: de_all(:,:)                  ! (nE,nsppol) differences between selected energies
260      real(dp),allocatable :: veloc_sq0(:,:,:)             ! (3,nsppol,nT)
261      real(dp),allocatable :: veloc_sq(:,:,:)              ! (3,nsppol,nE)
262 
263      real(dp),allocatable :: el_veloc(:,:,:,:)        ! nkpt nband 3 nsppol
264 ! the 9 = 3x3 is for the full tensorial transport coefficients
265      real(dp),allocatable :: gamma_qpt_tr(:,:,:,:,:)    ! 2 9 branches**2 nsppol qpt
266      real(dp),allocatable :: gamma_qpt_trin(:,:,:,:,:)  !idem
267      real(dp),allocatable :: gamma_qpt_trout(:,:,:,:,:) !idem
268 
269      real(dp),allocatable :: gamma_rpt_tr(:,:,:,:,:,:,:)    !idem
270      real(dp),allocatable :: gamma_rpt_trin(:,:,:,:,:)  !idem
271      real(dp),allocatable :: gamma_rpt_trout(:,:,:,:,:) !idem
272 
273      real(dp),allocatable :: a2f_1d_tr(:,:,:,:,:,:)           ! nfreq 9 nsppol 4 n_pair ntemp
274      real(dp),allocatable :: a2f_1d_trin(:,:,:)
275      real(dp),allocatable :: a2f_1d_trout(:,:,:)
276 
277      real(dp),allocatable :: FSelecveloc_sq(:,:)       ! 3 nsppol
278 
279      real(dp),allocatable :: tmp_gkk_intweight(:,:,:,:)
280      real(dp),allocatable :: tmp_gkk_intweight1(:,:,:)
281      real(dp),allocatable :: tmp_gkk_intweight2(:,:,:)
282 
283      real(dp),allocatable :: tmp_velocwtk(:,:,:,:,:)
284      real(dp),allocatable :: tmp_velocwtk1(:,:,:,:)
285      real(dp),allocatable :: tmp_velocwtk2(:,:,:,:)
286 
287      real(dp),allocatable :: tmp_vvelocwtk(:,:,:,:,:,:)
288      real(dp),allocatable :: tmp_vvelocwtk1(:,:,:,:,:)
289      real(dp),allocatable :: tmp_vvelocwtk2(:,:,:,:,:)
290 
291   end type elph_tr_type
292 
293  public :: elph_tr_ds_clean

defs_elphon/elph_type [ Types ]

[ Top ] [ defs_elphon ] [ Types ]

NAME

 elph_type

FUNCTION

 elph_type contains data and dimensions for the kpoints near the
 fermi surface and the $g_{k k+q}$ matrix elements

SOURCE

115   type,public :: elph_type
116 
117    type(elph_kgrid_type) :: k_phon               ! object for k-grid of phonon calculation
118    type(elph_kgrid_type) :: k_fine               ! object for fine k-grid for FS integration
119 
120    integer :: natom,nbranch,nFSband,nband
121    integer :: minFSband,maxFSband                !Index of lower and upper bands used for the FS integration
122 
123    integer :: ngkkband                           !Number of bands kept in final gkk matrix elements:
124                                                  !either 1 if sum is performed immediately
125                                                  !or = nband if all elements are kept based on flag ep_keepbands
126 
127    integer :: nenergy                            ! number of points in energy
128                                                  !   space for the electron band energies
129                                                  !   energies from
130                                                  !   ef-nenergy*delta_e to
131                                                  !   ef+nenergy*delta_e
132 
133    integer :: n_pair                             ! number of pairs considered
134 
135    integer :: nqpt_full                          !number of q in full BZ
136    integer :: nqptirred                          !number of irred q-points
137 
138 
139    integer :: unita2f,unit_gkk2,unit_gkk_rpt
140    integer :: unitgkq                            !units for file output
141 
142    integer :: gkqwrite
143    integer :: gkk2write
144    integer :: gkk_rptwrite
145 
146    integer :: ep_scalprod                        !flag to perform the scalar product
147    integer :: symgkq                             !flag to symmetrize gkq matrix elements
148    integer :: ep_keepbands                       !flag to sum over bands or not
149    integer :: ep_lova                            ! 1 for lova, and 0 for general
150    integer :: tuniformgrid                       !flag to expect uniform grid of q or not
151    integer :: prtbltztrp                         !flag to output BoltzTraP input files
152 
153    integer :: na2f                               !dimensions and increments for a2F function
154    integer :: nsppol                             ! number of spin polarization channels
155    integer :: nspinor                            ! number of spinorial components
156    integer :: telphint                           ! flag for integration over the FS with 0=tetrahedra 1=gaussians
157    integer :: ep_nspline                         ! scale factor for spline interpolation in RTA
158    integer :: ntemper                            ! number of temperature points
159    integer :: use_k_fine                         ! flag for using fine k-grids for eigenvalues and velocities. 0=no 1=yes
160    integer :: ep_int_gkk                         ! flag for interpolate gkk(1) or gamma (0)
161    integer :: ep_b_min                           ! first band taken into account in FS integration (if telphint==2)
162    integer :: ep_b_max                           ! last band taken into account in FS integration (if telphint==2)
163    integer :: kptrlatt(3,3)                      ! kpoint grid generating vectors, as in abinit
164    integer :: kptrlatt_fine(3,3)                 ! kpoint grid generating vectors, for fine grid used in FS integration
165 
166    real(dp) :: delta_e                           ! step in electronic energies, around Fermi level
167    real(dp) :: omega_min,omega_max
168    real(dp) :: a2fsmear,domega
169    real(dp) :: nelect                            ! number of electrons per unit cell, eventually with extra charges for carriers in semiconductors.
170    real(dp) :: occ_factor                        ! normalization for integrals over FS, for num of spins, spinors, etc...
171 
172    real(dp) :: mustar                            ! mustar parameter
173    real(dp) :: fermie                            ! Fermi energy (Ha), either comes from wfk file or from anaddb input file
174    real(dp) :: elphsmear                         ! smearing width for gaussian integration or buffer in energy for
175                                                  ! calculations with tetrahedra (telphint=0)
176    real(dp) :: tempermin                         ! minimum temperature at which resistivity etc are calculated (in K)
177    real(dp) :: temperinc                         ! interval temperature grid on which resistivity etc are calculated (in K)
178 
179    character(len=fnlen) :: elph_base_name        !base name for output files
180 
181    integer,allocatable :: qirredtofull(:)            !mapping between the qpoints found in the GGK file
182                                                  !and the array of qpoints generated by the code
183 
184    real(dp),allocatable :: wtq(:)                    !weight for each qpoint in the full grid spqt
185                                                  !if a point is not in the IBZ ==>  wtq=0
186                                                  !MG we can also use indqpt
187 
188    real(dp),allocatable :: n0(:)                     !DOS at the Fermi level (states/Ha/spin)
189    real(dp),allocatable :: qpt_full(:,:)             !special q points obtained by the Monkhorst & Pack method,
190                                                  !in reduced coordinates
191 
192 
193    real(dp),allocatable :: gkk_intweight(:,:,:)      ! (nFSband,nkpt_fine,nsppol)
194                                                  !integration weights for gkk matrix elements on FS:
195                                                  !if ep_keepbands == 0 all are 1
196                                                  !if ep_keepbands == 1 then = to wtk_phon in elphon
197                                                  !DOES NOT INCLUDE FACTOR OF 1/nkpt_phon
198 
199    real(dp),allocatable :: gkk_velocwtk(:,:,:)      ! (nFSband,nkpt_fine,nsppol)
200 
201    real(dp),allocatable :: gkk_vvelocwtk(:,:,:)      ! (nFSband,nkpt_fine,nsppol)
202 
203    real(dp),allocatable :: gkk_qpt(:,:,:,:,:,:)      ! (2, ngkkband*ngkkband, nbranch*nbranch, nkpt_phon, nsppol, nqptirred)
204                                                  !Now gkq contains gkk2 matrices on basic qpts,
205                                                  !summed over bands if ngkkband==1
206 
207 
208    real(dp),allocatable :: gkk_rpt(:,:,:,:,:,:)      ! (2, ngkkband**2, nbranch**2, nkpt_phon, nsppol, nrpt)
209                                                  !For the moment, gkk_rpt in memory is out of the question
210    real(dp),allocatable :: gkk2(:,:,:,:,:,:)         ! (nbranch, ngkkband,ngkkband, nkpt_phon, nkpt_phon, nsppol)
211 
212    real(dp),allocatable :: gamma_qpt(:,:,:,:)        !gamma matrices integrated over kpoint coeff
213                                                  !  and bands: still depends on qpt
214                                                  ! dims= 2, nbranch**2, nsppol, nqpt
215    real(dp),allocatable :: gamma_rpt(:,:,:,:)
216                                                  ! dims= 2, nbranch**2, nsppol, nrpt
217 !NOTE: choice to put nsppol before or after nqpt is a bit arbitrary
218 !   abinit uses nband,nkpt,nsppol, but here for convenience nkpt_phon,nsppol,nqpt
219 !   as interpolation is on qpt
220 
221    real(dp),allocatable :: phfrq(:,:)                !phonon frequencies
222    real(dp),allocatable :: a2f(:,:,:)                !a2f function
223 
224    real(dp),allocatable :: qgrid_data(:,:,:,:)       !e-ph values calculated over the irreducible part of the q-grid:
225                                                  !first entry  =  index of the q-point,
226                                                  !second index =  branch index
227                                                  !the third slice contains the frequency, the linewidth and lambda(q,nu)
228                                                  !for that particular phonon mode
229                                                  ! dims= nqptirred,elph_ds%nbranch,nsppol,3
230 
231  end type elph_type
232 
233  public :: elph_ds_clean

defs_elphon/gam_mult_displ [ Functions ]

[ Top ] [ defs_elphon ] [ Functions ]

NAME

 gam_mult_displ

FUNCTION

 This routine takes the bare gamma matrices and multiplies them
  by the displ_red matrices (related to the scalprod variable)

INPUTS

   nbranch = number of phonon branches (3*natom)
   displ_red = phonon mode displacement vectors in reduced coordinates.
   gam_bare = bare gamma matrices before multiplication

OUTPUT

   gam_now = output gamma matrices multiplied by displacement matrices

PARENTS

      get_tau_k,m_phgamma,mka2f,mka2f_tr,mka2f_tr_lova,mkph_linwid,nmsq_gam
      nmsq_gam_sumfs,nmsq_pure_gkk,nmsq_pure_gkk_sumfs,normsq_gkq

CHILDREN

      zgemm

SOURCE

792 subroutine gam_mult_displ(nbranch, displ_red, gam_bare, gam_now)
793 
794 
795 !This section has been created automatically by the script Abilint (TD).
796 !Do not modify the following lines by hand.
797 #undef ABI_FUNC
798 #define ABI_FUNC 'gam_mult_displ'
799 !End of the abilint section
800 
801  implicit none
802 
803 !Arguments -------------------------------
804  integer, intent(in)  :: nbranch
805  real(dp), intent(in)  :: displ_red(2,nbranch,nbranch)
806  real(dp), intent(in)  :: gam_bare(2,nbranch,nbranch)
807  real(dp), intent(out) :: gam_now(2,nbranch,nbranch)
808 
809 !Local variables -------------------------
810  real(dp) :: zgemm_tmp_mat(2,nbranch,nbranch)
811 
812 ! *********************************************************************
813 
814  gam_now = zero
815 
816  call zgemm('c','n',nbranch,nbranch,nbranch,cone,&
817 & displ_red,nbranch,gam_bare,&
818 & nbranch,czero,zgemm_tmp_mat,nbranch)
819 
820  call zgemm('n','n',nbranch,nbranch,nbranch,cone,&
821 & zgemm_tmp_mat,nbranch,displ_red,&
822 & nbranch,czero,gam_now,nbranch)
823 
824 end subroutine gam_mult_displ