TABLE OF CONTENTS


ABINIT/m_paw_dmft [ Modules ]

[ Top ] [ Modules ]

NAME

  m_paw_dmft

FUNCTION

COPYRIGHT

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

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

24 #if defined HAVE_CONFIG_H
25 #include "config.h"
26 #endif
27 
28 #include "abi_common.h"
29 
30 MODULE m_paw_dmft
31 
32  use defs_basis
33  use defs_datatypes
34  use defs_abitypes
35  use m_CtqmcInterface
36  use m_errors
37  use m_abicore
38  use m_xmpi
39  use m_data4entropyDMFT
40 
41  use m_io_tools,  only : open_file
42  use m_pawtab,    only : pawtab_type
43 
44  implicit none
45 
46  private
47 
48  public :: init_dmft
49  public :: init_sc_dmft
50  public :: construct_nwli_dmft
51  public :: destroy_dmft
52  public :: destroy_sc_dmft
53  public :: print_dmft
54  public :: print_sc_dmft
55  public :: saveocc_dmft
56  public :: readocc_dmft

m_paw_dmft/construct_nwli_dmft [ Functions ]

[ Top ] [ m_paw_dmft ] [ Functions ]

NAME

 construct_nwli_dmft

FUNCTION

  Compute linear frequencies

INPUTS

  paw_dmft=structure for dmft
  nwli=number of linear frequencies

 OUTPUTS
  omegali(1:nwli)=computed frequencies

PARENTS

      m_green,m_paw_dmft

CHILDREN

      wrtout

SOURCE

840 subroutine construct_nwli_dmft(paw_dmft,nwli,omega_li)
841 
842 
843 !This section has been created automatically by the script Abilint (TD).
844 !Do not modify the following lines by hand.
845 #undef ABI_FUNC
846 #define ABI_FUNC 'construct_nwli_dmft'
847 !End of the abilint section
848 
849  implicit none
850 
851  type(paw_dmft_type), intent(in) :: paw_dmft
852  integer, intent(in) :: nwli
853  real(dp), intent(out) :: omega_li(:)
854  !fortran2003 ?
855  !real(dp), allocatable, intent(inout) :: omega_li(:)
856  integer :: ifreq
857  real(dp) :: factor
858  character(len=100) :: message
859 
860 ! if (allocated(omega_li)) then
861    if (size(omega_li) .ne. nwli) then
862      write(message,'(2a,i8,a,i8)') ch10, "Number of linear frequencies asked is", &
863        &    nwli, "whereas dimension of array omega_li is", size(omega_li)
864      MSG_BUG(message)
865 !     ABI_DEALLOCATE(omega_li)
866 !     ABI_ALLOCATE(omega_li,(nwli))
867 !     write(*,*) "RESIZE"
868 !     call flush(6)
869    endif
870 !     write(*,*) "NOTHING"
871 !     call flush(6)
872 ! else
873 !     write(*,*) "ALLOCATE"
874 !     call flush(6)
875 !   ABI_ALLOCATE(omega_li,(nwli))
876 ! endif
877 
878 ! Set up linear frequencies
879  factor = pi*paw_dmft%temp
880  do ifreq=1,nwli
881    omega_li(ifreq)=factor*(real(2*ifreq-1,kind=dp))
882    ! (2(ifreq-1)+1 = 2ifreq-1
883  enddo
884 end subroutine construct_nwli_dmft

m_paw_dmft/construct_nwlo_dmft [ Functions ]

[ Top ] [ m_paw_dmft ] [ Functions ]

NAME

 construct_nwlo_dmft

FUNCTION

  Allocate log frequencies if used and compute them as well as their weight

INPUTS

  paw_dmft=structure for dmft calculation

PARENTS

      m_paw_dmft

CHILDREN

      wrtout

SOURCE

 906 !! NOTE
 907 !! The part of the code which deals
 908 !! with the use of logarithmic frequencies
 909 !! is a modification of the GNU GPL
 910 !! code available on http://dmft.rutgers.edu/ and
 911 !! described in the  RMP paper written by
 912 !! G.Kotliar,  S.Y.Savrasov, K.Haule, V.S.Oudovenko, O.Parcollet, C.A.Marianetti.
 913 !!
 914 
 915 subroutine construct_nwlo_dmft(paw_dmft)
 916  use m_splines
 917 
 918 !This section has been created automatically by the script Abilint (TD).
 919 !Do not modify the following lines by hand.
 920 #undef ABI_FUNC
 921 #define ABI_FUNC 'construct_nwlo_dmft'
 922 !End of the abilint section
 923 
 924  implicit none
 925 
 926  type(paw_dmft_type), intent(inout) :: paw_dmft
 927  integer :: cubic_freq
 928  integer :: ifreq,ifreq2
 929  ! for parallel
 930  integer :: myproc, nproc, spacecomm
 931  integer :: deltaw, residu, omegaBegin, omegaEnd
 932  ! end
 933  character(len=500) :: message
 934  real(dp) :: deltaomega,expfac,omegamaxmin,prefacexp,AA,BB,CC,nlin,nlog,t1
 935  real(dp) :: wl
 936  complex(dpc):: ybcbeg,ybcend
 937  integer, allocatable :: select_log(:)
 938  real(dp), allocatable :: omega_lo_tmp(:)
 939  real(dp), allocatable :: omega_li(:)
 940  real(dp), allocatable :: wgt_wlo(:)
 941  complex(dpc), allocatable :: tospline_lo(:), splined_li(:),ysplin2_lo(:)
 942 
 943  ABI_ALLOCATE(omega_lo_tmp,(paw_dmft%dmft_nwlo))
 944  ABI_ALLOCATE(wgt_wlo,(paw_dmft%dmft_nwlo))
 945 
 946 !==  Variables for DMFT related to frequencies
 947 ! the part of the code which deals
 948 ! with the use of logarithmic frequencies
 949 ! is a modification of the GNU GPL
 950 ! code available on http://dmft.rutgers.edu/ and
 951 ! described in the  RMP paper written by
 952 ! G.Kotliar, S.Y.Savrasov, K.Haule, V.S.Oudovenko, O.Parcollet, C.A.Marianetti
 953 
 954 !========================================
 955 !== construct log. freq.
 956  if(paw_dmft%dmft_log_freq==1) then
 957 !=======================================
 958    cubic_freq=0
 959    !omegamaxmin=paw_dmft%omega_li(paw_dmft%dmft_nwli)-paw_dmft%omega_li(paw_dmft%dmftqmc_l+1)
 960    omegamaxmin=pi*paw_dmft%temp*two*real(paw_dmft%dmft_nwli-paw_dmft%dmftqmc_l-1,kind=dp)
 961 
 962    if(cubic_freq==1) then
 963 
 964      if (paw_dmft%dmft_solv .eq. 5 ) then
 965        write(message, '(2a)') ch10, "Warning : Cubish Mesh not tested with CT-QMC"
 966        MSG_WARNING(message)
 967      end if
 968 !  ------------  CUBIC MESH MESH
 969 !    useless
 970      nlin=dble(paw_dmft%dmft_nwli)
 971      nlog=dble(paw_dmft%dmft_nwlo)
 972      AA=((nlin-one)/nlin/(nlog**2-one)-one/(three*nlin))/((nlog**3-one)/(nlog**2-one)-seven/three)
 973      BB=(one/nlin - seven*AA)/three
 974      CC=-AA-BB
 975 !    AA=((nlin-one)/nlin/(nlog-one)-one/(nlin))/((nlog**2-one)/(nlog-one)-three)
 976 !    BB=(one/nlin - three*AA)
 977 !    CC=-AA-BB
 978      write(message, '(a,16x,2(2x,a))') ch10,"  Cubic Mesh Parameters are"
 979      call wrtout(std_out,message,'COLL')
 980      write(message, '(3x,a,3(2x,e13.5))') "AA,BB,CC",AA,BB,CC
 981      call wrtout(std_out,message,'COLL')
 982      do ifreq=1,paw_dmft%dmft_nwlo
 983        t1=dble(ifreq)
 984        !omega_lo_tmp(ifreq)=(AA*t1**3+BB*t1**2+CC)*omegamaxmin+paw_dmft%omega_li(1)
 985        omega_lo_tmp(ifreq)=(AA*t1**3+BB*t1**2+CC)*omegamaxmin+paw_dmft%temp*pi
 986 !       paw_dmft%omega_lo(ifreq)=(AA*t1**2+BB*t1+CC)*omegamaxmin+paw_dmft%omega_li(1)
 987 !     write(69,*) paw_dmft%omega_lo(ifreq),0.5
 988      enddo
 989    else
 990      if(paw_dmft%dmft_solv<4) paw_dmft%dmftqmc_l=0
 991 
 992 !   ------------  LOGARITHMIC MESH
 993      deltaomega=0.5_dp
 994      expfac=log(omegamaxmin/deltaomega)/(float(paw_dmft%dmft_nwlo-paw_dmft%dmftqmc_l-1)/two)
 995      prefacexp=omegamaxmin/(exp(expfac*float(paw_dmft%dmft_nwlo-paw_dmft%dmftqmc_l-1))-one)
 996      ABI_ALLOCATE(select_log,(paw_dmft%dmft_nwlo))
 997      select_log=0
 998 
 999 !   ------------ IMPOSE LINEAR MESH for w < 2*w_n=(2*l-1)pi/beta
1000 !         Check variables (Already done in chkinp if dmft_solv==5)
1001      if (paw_dmft%dmftqmc_l .gt. paw_dmft%dmft_nwlo) then
1002        write(message, '(a,a,i6)' )ch10,&
1003 &       ' ERROR: dmft_nwlo has to be at leat equal to 2xdmftqmc_l :',2*paw_dmft%dmftqmc_l
1004        MSG_ERROR(message)
1005      end if
1006 !         End Check
1007 
1008      call construct_nwli_dmft(paw_dmft,paw_dmft%dmftqmc_l,omega_lo_tmp(1:paw_dmft%dmftqmc_l))
1009      select_log(1:paw_dmft%dmftqmc_l) = (/ (ifreq,ifreq=1,paw_dmft%dmftqmc_l) /)
1010 
1011      !do ifreq=1,paw_dmft%dmftqmc_l
1012      !  omega_lo_tmp(ifreq)=(two*DBLE(ifreq-1)+one)*pi*paw_dmft%temp
1013      !  select_log(ifreq)=ifreq
1014      !enddo
1015 
1016 !   ------------ COMPLETE FREQUENCIES WITH LOG MESH
1017      wl = paw_dmft%temp*pi*real(2*paw_dmft%dmftqmc_l+1,kind=dp)
1018      do ifreq=1,paw_dmft%dmft_nwlo-paw_dmft%dmftqmc_l
1019        !omega_lo_tmp(ifreq+paw_dmft%dmftqmc_l)=prefacexp*(exp(expfac*float(ifreq-1))-one)+paw_dmft%omega_li(paw_dmft%dmftqmc_l+1)
1020        omega_lo_tmp(ifreq+paw_dmft%dmftqmc_l)=prefacexp*(exp(expfac*real(ifreq-1,kind=dp))-one) + wl
1021 !       -------- Impose that the each frequency of the logarithmic mesh is on a Matsubara frequency
1022 ! FIXME : This may be done for all solver, not only for QMCs
1023        if(paw_dmft%dmft_solv>=4) then
1024          ! Compute the integer "n" of iwn
1025          ifreq2 = nint((omega_lo_tmp(ifreq+paw_dmft%dmftqmc_l)/(paw_dmft%temp*pi)-one)*half)
1026          ! compute freq
1027          omega_lo_tmp(ifreq+paw_dmft%dmftqmc_l)= (dble(ifreq2)*two+one)*pi*paw_dmft%temp
1028 
1029          if((ifreq2+1)>paw_dmft%dmft_nwli) then
1030            write(message, '(a,a,i8)' )ch10,&
1031 &          ' BUG: init_dmft,   dimension  of array select_log is about to be overflown',&
1032 &          (ifreq2+1)
1033            MSG_BUG(message)
1034          endif
1035          select_log(paw_dmft%dmftqmc_l+ifreq)=ifreq2+1
1036        endif
1037      enddo
1038 
1039 !       -------- Suppress duplicate frequencies
1040 ! FIXME : So this also should be done for all solver and remove useless
1041 ! frequencies
1042      if(paw_dmft%dmft_solv>=4) then
1043        ifreq2=1
1044        do ifreq=2,paw_dmft%dmft_nwlo-1
1045          if(select_log(ifreq2).ne.select_log(ifreq)) then
1046            ifreq2=ifreq2+1
1047            omega_lo_tmp(ifreq2)=omega_lo_tmp(ifreq)
1048            select_log(ifreq2)=select_log(ifreq)
1049          endif
1050        enddo
1051        paw_dmft%dmft_nwlo=ifreq2+1
1052      endif
1053 
1054      omega_lo_tmp(1)=paw_dmft%temp*pi
1055      omega_lo_tmp(paw_dmft%dmft_nwlo)=paw_dmft%temp*pi*real(2*paw_dmft%dmft_nwli-1,kind=dp)
1056   endif
1057 
1058 !=======================
1059 !== construct weight for log. freq.
1060 !=======================
1061   ABI_ALLOCATE(tospline_lo,(paw_dmft%dmft_nwlo))
1062   ABI_ALLOCATE(splined_li,(paw_dmft%dmft_nwli))
1063   ABI_ALLOCATE(ysplin2_lo,(paw_dmft%dmft_nwlo))
1064   if (allocated(omega_li)) then
1065     ABI_DEALLOCATE(omega_li)
1066   endif
1067   ABI_ALLOCATE(omega_li,(1:paw_dmft%dmft_nwli))
1068   call construct_nwli_dmft(paw_dmft,paw_dmft%dmft_nwli,omega_li)
1069 
1070   !Parallelisation over frequencies!
1071   ! ============= Set up =============
1072   myproc = paw_dmft%myproc
1073   nproc = paw_dmft%nproc
1074   spacecomm = paw_dmft%spacecomm
1075   deltaw = paw_dmft%dmft_nwlo / nproc
1076   residu = paw_dmft%dmft_nwlo - nproc*deltaw
1077   if ( myproc .LT. nproc - residu ) then
1078     omegaBegin = 1 + myproc*deltaw
1079     omegaEnd   = (myproc + 1)*deltaw
1080   else
1081     omegaBegin = 1 + myproc*(deltaw + 1) -nproc + residu
1082     omegaEnd = omegaBegin + deltaw
1083   end if
1084   wgt_wlo(:)=zero
1085   ! ============= END Set up =============
1086   do ifreq=omegaBegin,omegaEnd
1087     tospline_lo=cmplx(0_dp,0_dp,kind=dp)
1088 !    do ifreq1=1,paw_dmft%dmft_nwlo
1089     tospline_lo(ifreq)=cmplx(1_dp,0_dp,kind=dp)
1090 !    tospline_lo(ifreq1)=ifreq1**2-ifreq1
1091 !    enddo
1092     splined_li=cmplx(0_dp,0_dp,kind=dp)
1093 !    ybcbeg=cmplx(one/tol16**2,zero)
1094 !    ybcend=cmplx(one/tol16**2,zero)
1095     ybcbeg=czero
1096     ybcend=czero
1097 
1098 
1099 !==         spline delta function
1100     call spline_complex( omega_lo_tmp, tospline_lo, paw_dmft%dmft_nwlo, &
1101    & ybcbeg, ybcend, ysplin2_lo)
1102 !   do ifreq1=1,paw_dmft%dmft_nwlo
1103 !    write(6588,*) paw_dmft%omega_lo(ifreq1),ysplin2_lo(ifreq1)
1104 !   enddo
1105 
1106     call splint_complex( paw_dmft%dmft_nwlo, omega_lo_tmp, tospline_lo,&
1107    & ysplin2_lo, paw_dmft%dmft_nwli, omega_li, splined_li)
1108 
1109 !==         accumulate weights
1110     wgt_wlo(ifreq)=zero
1111     do ifreq2=1,paw_dmft%dmft_nwli
1112       wgt_wlo(ifreq)=wgt_wlo(ifreq)+real(splined_li(ifreq2),kind=dp)
1113     enddo
1114 ! do ifreq1=1,paw_dmft%dmft_nwlo
1115 !  write(6688,*) paw_dmft%omega_lo(ifreq1),tospline_lo(ifreq1)
1116 ! enddo
1117 ! do ifreq1=1,paw_dmft%dmft_nwli
1118 !  write(6788,*) paw_dmft%omega_li(ifreq1),splined_li(ifreq1)
1119   enddo
1120   ! ============= Gatherall  =============
1121   call xmpi_sum(wgt_wlo, spacecomm, residu)
1122   ! ============= END Gatherall ==========
1123   ! end parallelisation over frequencies
1124 
1125   ABI_DEALLOCATE(tospline_lo)
1126   ABI_DEALLOCATE(splined_li)
1127   ABI_DEALLOCATE(ysplin2_lo)
1128 ! if(abs(dtset%pawprtvol)>=3) then
1129     write(message, '(a,18x,2(2x,a21))') ch10,"Log. Freq","weight"
1130     call wrtout(std_out,message,'COLL')
1131     do ifreq=1,paw_dmft%dmft_nwlo
1132       write(message, '(3x,a9,i6,2(2x,e21.14))') "--ifreq--",ifreq,omega_lo_tmp(ifreq),wgt_wlo(ifreq)
1133       call wrtout(std_out,message,'COLL')
1134     enddo
1135     write(message, '(3x,a,i6)') "  Total number of log frequencies is", paw_dmft%dmft_nwlo
1136     call wrtout(std_out,message,'COLL')
1137     ifreq2 = 1
1138     do ifreq=1,min(30,paw_dmft%dmft_nwlo)
1139       write(message, '(3x,a9,i6,2(2x,e21.14))') "--ifreq--",ifreq,omega_li(ifreq)
1140       call wrtout(std_out,message,'COLL')
1141       if ( select_log(ifreq2) .eq. ifreq ) then
1142         write(message, '(3x,a,i4,2(2x,i5))') "--sel_log",1
1143         ifreq2 = ifreq+1
1144       else
1145         write(message, '(3x,a,i4,2(2x,i5))') "--sel_log",0
1146       end if
1147       call wrtout(std_out,message,'COLL')
1148     enddo
1149     write(message, '(3x,2a)') "--ifreq--","..."
1150     call wrtout(std_out,message,'COLL')
1151     write(message, '(3x,a,i6,2(2x,e13.5))') "--ifreq--",paw_dmft%dmft_nwli,omega_li(paw_dmft%dmft_nwli)
1152     call wrtout(std_out,message,'COLL')
1153 !   endif
1154    ABI_DEALLOCATE(select_log)
1155    ABI_DEALLOCATE(omega_li)
1156 
1157 !=========================================================
1158 !== do not construct log. freq. and use linear frequencies
1159  else
1160 !=========================================================
1161    write(message, '(a,10x,2(2x,a))') ch10,"   Use of linear frequencies for DMFT calculation"
1162    call wrtout(std_out,message,'COLL')
1163    call construct_nwli_dmft(paw_dmft,paw_dmft%dmft_nwli,omega_lo_tmp)
1164    wgt_wlo=one
1165  endif
1166 
1167  ! Should be check but since type definition does not initialize pointer with
1168  ! =>null() (fortran95 and later) it produces conditional jump in valgrind
1169  !if ( associated(paw_dmft%omega_lo) ) then
1170  !  ABI_DEALLOCATE(paw_dmft%omega_lo)
1171  !endif
1172  !if ( associated(paw_dmft%wgt_wlo) ) then
1173  !  ABI_DEALLOCATE(paw_dmft%wgt_wlo)
1174  !endif
1175  ABI_ALLOCATE(paw_dmft%omega_lo,(paw_dmft%dmft_nwlo))
1176  ABI_ALLOCATE(paw_dmft%wgt_wlo,(paw_dmft%dmft_nwlo))
1177  paw_dmft%omega_lo(1:paw_dmft%dmft_nwlo) = omega_lo_tmp(1:paw_dmft%dmft_nwlo)
1178  paw_dmft%wgt_wlo(1:paw_dmft%dmft_nwlo) = wgt_wlo(1:paw_dmft%dmft_nwlo)
1179  ABI_DEALLOCATE(omega_lo_tmp)
1180  ABI_DEALLOCATE(wgt_wlo)
1181 end subroutine construct_nwlo_dmft

m_paw_dmft/destroy_dmft [ Functions ]

[ Top ] [ m_paw_dmft ] [ Functions ]

NAME

 destroy_dmft

FUNCTION

  deallocate some variables related to paw_dmft

INPUTS

  paw_dmft

OUTPUT

PARENTS

      outscfcv,vtorho

CHILDREN

      wrtout

SOURCE

1204 subroutine destroy_dmft(paw_dmft)
1205 
1206 
1207 !This section has been created automatically by the script Abilint (TD).
1208 !Do not modify the following lines by hand.
1209 #undef ABI_FUNC
1210 #define ABI_FUNC 'destroy_dmft'
1211 !End of the abilint section
1212 
1213  implicit none
1214 
1215 !Arguments ------------------------------------
1216 !scalars
1217  type(paw_dmft_type),intent(inout) :: paw_dmft
1218 
1219 !Local variables-------------------------------
1220  integer :: iatom
1221 
1222 ! *********************************************************************
1223 
1224    if (paw_dmft%dmft_solv == 5 .and. allocated(paw_dmft%hybrid)) then
1225      do iatom=1, size(paw_dmft%hybrid) !paw_dmft%natom
1226        !if(paw_dmft%lpawu(iatom)/=-1) then
1227          call ctqmcinterface_finalize(paw_dmft%hybrid(iatom))
1228        !endif
1229      enddo
1230      ABI_DATATYPE_DEALLOCATE(paw_dmft%hybrid)
1231    endif
1232    if (allocated(paw_dmft%psichi))  then
1233      ABI_DEALLOCATE(paw_dmft%psichi)
1234    end if
1235 !   paw_dmft%wtk is only an explicit pointer =>dtset%wtk
1236 !   if (associated(paw_dmft%wtk)) deallocate(paw_dmft%wtk)
1237    paw_dmft%wtk => null()
1238    paw_dmft%fixed_self => null()
1239    if (allocated(paw_dmft%eigen_lda))  then
1240      ABI_DEALLOCATE(paw_dmft%eigen_lda)
1241    endif
1242    if (associated(paw_dmft%omega_lo))  then
1243      ABI_DEALLOCATE(paw_dmft%omega_lo)
1244    end if
1245    if (associated(paw_dmft%omega_r))  then
1246      ABI_DEALLOCATE(paw_dmft%omega_r)
1247    end if
1248    if (associated(paw_dmft%wgt_wlo))  then
1249      ABI_DEALLOCATE(paw_dmft%wgt_wlo)
1250    end if
1251    if (allocated(paw_dmft%lpawu))  then
1252      ABI_DEALLOCATE(paw_dmft%lpawu)
1253    end if
1254 
1255 end subroutine destroy_dmft

m_paw_dmft/destroy_sc_dmft [ Functions ]

[ Top ] [ m_paw_dmft ] [ Functions ]

NAME

 destroy_sc_dmft

FUNCTION

  deallocate paw_dmft

INPUTS

  paw_dmft

OUTPUT

PARENTS

      gstate

CHILDREN

      wrtout

SOURCE

1278 subroutine destroy_sc_dmft(paw_dmft)
1279 
1280 
1281 !This section has been created automatically by the script Abilint (TD).
1282 !Do not modify the following lines by hand.
1283 #undef ABI_FUNC
1284 #define ABI_FUNC 'destroy_sc_dmft'
1285 !End of the abilint section
1286 
1287  implicit none
1288 
1289 !Arguments ------------------------------------
1290 !scalars
1291  type(paw_dmft_type),intent(inout) :: paw_dmft
1292 
1293 !Local variables-------------------------------
1294  character(len=500) :: message
1295 
1296 ! *********************************************************************
1297 
1298  if (( .not. allocated(paw_dmft%occnd) .or. .not. allocated(paw_dmft%band_in) &
1299 &  .or. .not. allocated(paw_dmft%include_bands) .or. .not. allocated(paw_dmft%exclude_bands)) &
1300 &  .and. paw_dmft%use_dmft == 1 )  then
1301   write(message, '(a,a,a)' )&
1302 &  '  an array is not allocated and is not deallocated with use_dmft==1 ',ch10, &
1303 &  '  Action : check the code'
1304   MSG_WARNING(message)
1305  endif
1306  if ( allocated(paw_dmft%occnd) )          then
1307    ABI_DEALLOCATE(paw_dmft%occnd)
1308  end if
1309  if ( allocated(paw_dmft%band_in) )        then
1310    ABI_DEALLOCATE(paw_dmft%band_in)
1311  end if
1312  if ( allocated(paw_dmft%include_bands) )  then
1313    ABI_DEALLOCATE(paw_dmft%include_bands)
1314  end if
1315  if ( allocated(paw_dmft%exclude_bands) )  then
1316    ABI_DEALLOCATE(paw_dmft%exclude_bands)
1317  end if
1318 
1319 
1320 end subroutine destroy_sc_dmft

m_paw_dmft/init_dmft [ Functions ]

[ Top ] [ m_paw_dmft ] [ Functions ]

NAME

 init_dmft

FUNCTION

  Allocate variables and setup lda hamiltonian and related data
  (init_sc_dmft has to been called before)

INPUTS

  dmatpawu   = fixed occupation matrix of correlated orbitals
  eigen      = LDA eigenvalues
  fermie_lda = LDA Fermi level
  psichi     = <chi|Psi> projection of KS states over atomic !wavefunction
  nkpt       = number of k-points
  nsppol     = number of spin polarisation
  nspinor    = number of spinorial component

PARENTS

      outscfcv,vtorho

CHILDREN

      wrtout

SOURCE

528 !! NOTE
529 !! The part of the code which deals
530 !! with the use of logarithmic frequencies
531 !! is a modification of the GNU GPL
532 !! code available on http://dmft.rutgers.edu/ and
533 !! described in the  RMP paper written by
534 !! G.Kotliar,  S.Y.Savrasov, K.Haule, V.S.Oudovenko, O.Parcollet, C.A.Marianetti.
535 !!
536 
537 subroutine init_dmft(dmatpawu, dtset, fermie_lda, fnametmp_app, nspinor, paw_dmft, pawtab, psps, typat)
538 
539  use defs_abitypes
540  use m_splines
541  !use m_CtqmcInterface
542 
543 !This section has been created automatically by the script Abilint (TD).
544 !Do not modify the following lines by hand.
545 #undef ABI_FUNC
546 #define ABI_FUNC 'init_dmft'
547 !End of the abilint section
548 
549  implicit none
550 
551 !Arguments ------------------------------------
552 !scalars
553  integer, intent(in)  :: nspinor
554  real(dp), intent(in) :: fermie_lda
555 !type
556  type(pseudopotential_type), intent(in) :: psps
557  type(dataset_type),target,intent(in) :: dtset
558  type(pawtab_type),intent(in)  :: pawtab(psps%ntypat*psps%usepaw)
559  type(paw_dmft_type),intent(inout) :: paw_dmft
560  character(len=fnlen), intent(in) :: fnametmp_app
561 !arrays
562  integer,intent(in) :: typat(dtset%natom)
563  real(dp),intent(in),target :: dmatpawu(:,:,:,:)
564 !Local variables ------------------------------------
565  integer :: ikpt,isym,itypat,nsppol,mbandc,maxlpawu, iatom, ifreq
566  integer :: nflavor
567  real(dp) :: sumwtk
568  character(len=500) :: message
569  real(dp) :: unit_e,step
570 ! *********************************************************************
571 
572  nsppol = dtset%nsppol
573  if(dtset%ucrpa==0) then
574  write(message,'(6a)') ch10,' ====================================', &
575 &                      ch10,' =====  Start of DMFT calculation', &
576 &                      ch10,' ===================================='
577  else if(dtset%ucrpa>0) then
578  write(message,'(6a)') ch10,' ============================================================', &
579 &                      ch10,' =====  Initialize construction of Wannier in DMFT routines',&
580 &                      ch10,' ============================================================'
581  endif
582  call wrtout(std_out,message,'COLL')
583 
584  unit_e=2_dp
585 !=======================
586 !==  Check sym
587 !=======================
588  do isym=1,dtset%nsym
589    if(dtset%symafm(isym)<0) then
590      message = 'symafm negative is not implemented in DMFT '
591      MSG_ERROR(message)
592    endif
593  enddo
594 
595 ! paw_dmft%dmft_mag=0
596 ! do iatom=1,dtset%natom
597 !   do  ii=1,3
598 !     if ( dtset(ii,iatom) > 0.001 ) paw_dmft%dmft_mag=1
599 !   enddo
600 ! enddo
601 
602 !=======================
603 !==  Define integers and reals
604 !=======================
605  paw_dmft%fermie_lda=fermie_lda ! in Ha
606  paw_dmft%fermie= fermie_lda
607  if(nspinor==1) then
608    if(paw_dmft%nsppol==2) then
609      paw_dmft%nelectval= dtset%nelect-float(paw_dmft%dmftbandi-1)*paw_dmft%nsppol
610    else if(paw_dmft%nsppol==1) then
611      paw_dmft%nelectval= dtset%nelect-float(paw_dmft%dmftbandi-1)*paw_dmft%nsppol*2
612    endif
613  else if (nspinor==2) then
614    paw_dmft%nelectval= dtset%nelect-float(paw_dmft%dmftbandi-1)*paw_dmft%nsppol
615  endif
616  paw_dmft%filapp= fnametmp_app
617  paw_dmft%natpawu=dtset%natpawu
618  paw_dmft%natom=dtset%natom
619  paw_dmft%temp=dtset%tsmear!*unit_e
620  paw_dmft%dmft_iter=dtset%dmft_iter
621  paw_dmft%dmft_dc=dtset%dmft_dc
622  !paw_dmft%idmftloop=0
623  paw_dmft%prtvol = dtset%prtvol
624  paw_dmft%prtdos = dtset%prtdos
625  paw_dmft%dmft_tolfreq = dtset%dmft_tolfreq
626  paw_dmft%dmft_lcpr = dtset%dmft_tollc
627 
628 !=======================
629 !==  Fixed self for input
630 !=======================
631  paw_dmft%use_fixed_self=dtset%usedmatpu
632  paw_dmft%fixed_self=>dmatpawu
633 
634 !=======================
635 !==  Choose solver
636 !=======================
637  paw_dmft%dmft_solv=dtset%dmft_solv
638 !  0: LDA, no solver
639 !  1: LDA+U
640 ! -1: LDA+U but LDA values are not renormalized !
641 ! if((paw_dmft%dmft_solv==0.and.paw_dmft%prtvol>4).or.&
642 !&   (paw_dmft%dmft_solv>=-1.and.paw_dmft%dmft_solv<=2)) then
643 !   call wrtout(std_out,message,'COLL')
644 !   call wrtout(ab_out,message,'COLL')
645 ! endif
646 
647  if(paw_dmft%dmft_solv==0) then
648    do itypat=1,psps%ntypat
649      if(pawtab(itypat)%lpawu/=-1) then
650        if((pawtab(itypat)%upawu)>tol5.or.(pawtab(itypat)%jpawu)>tol5) then
651           write(message, '(2a,i5,2a,2e15.6)' )ch10,&
652 &          ' option dmft_solv=0 requires upaw=jpaw=0 for species',itypat,ch10,&
653 &          ' Value of upawu and jpawu are here',pawtab(itypat)%upawu,pawtab(itypat)%jpawu
654           MSG_ERROR(message)
655         endif
656      endif
657    enddo
658  endif
659 
660 ! todo_ab: why upaw and jpawu are not zero (on bigmac) if lpawu==-1 ?
661 ! if(paw_dmft%dmft_solv==0.and.&
662 !& (maxval(abs(pawtab(:)%upawu))>tol5.or.maxval(abs(pawtab(:)%jpawu))>tol5)) then
663 !   write(message, '(a,a,2f12.3)' )ch10,&
664 !&   ' option dmft_solv=0 requires upaw=jpaw=0',maxval(abs(pawtab(:)%upawu)),maxval(abs(pawtab(:)%jpawu))
665 !    MSG_WARNING(message)
666 ! endif
667 
668  paw_dmft%dmftcheck=dtset%dmftcheck
669 
670  if(paw_dmft%dmftcheck==-1) then
671    message = ' init_dmft: dmftcheck=-1 should not happend here'
672    MSG_BUG(message)
673  endif
674  paw_dmft%dmft_log_freq=1 ! use logarithmic frequencies.
675  if(paw_dmft%dmft_solv==6.or.paw_dmft%dmft_solv==7) then
676    paw_dmft%dmft_log_freq=0 ! do not use logarithmic frequencies.
677  endif
678  paw_dmft%dmft_nwli=dtset%dmft_nwli
679  if(paw_dmft%dmft_log_freq==1) then
680    paw_dmft%dmft_nwlo=dtset%dmft_nwlo
681  else
682    paw_dmft%dmft_nwlo=dtset%dmft_nwli
683  endif
684  paw_dmft%dmft_nwr=32
685  paw_dmft%dmft_rslf=dtset%dmft_rslf
686  paw_dmft%dmft_mxsf=dtset%dmft_mxsf
687  paw_dmft%dmftqmc_l=dtset%dmftqmc_l
688  paw_dmft%dmftqmc_n=dtset%dmftqmc_n
689  paw_dmft%dmftqmc_seed=dtset%dmftqmc_seed
690  paw_dmft%dmftqmc_therm=dtset%dmftqmc_therm
691  paw_dmft%dmftqmc_t2g=dtset%dmft_t2g
692  paw_dmft%dmftctqmc_basis =dtset%dmftctqmc_basis
693  paw_dmft%dmftctqmc_check =dtset%dmftctqmc_check
694  paw_dmft%dmftctqmc_correl=dtset%dmftctqmc_correl
695  paw_dmft%dmftctqmc_gmove =dtset%dmftctqmc_gmove
696  paw_dmft%dmftctqmc_grnns =dtset%dmftctqmc_grnns
697  paw_dmft%dmftctqmc_meas  =dtset%dmftctqmc_meas
698  paw_dmft%dmftctqmc_mrka  =dtset%dmftctqmc_mrka
699  paw_dmft%dmftctqmc_mov   =dtset%dmftctqmc_mov
700  paw_dmft%dmftctqmc_order =dtset%dmftctqmc_order
701  paw_dmft%dmftctqmc_triqs_nleg =dtset%dmftctqmc_triqs_nleg
702 
703  if ( paw_dmft%dmft_solv >= 4 ) then
704  write(message, '(a,a,i6)' )ch10,&
705 &   '=> Seed for QMC inside DMFT is dmftqmc_seed=',paw_dmft%dmftqmc_seed
706    call wrtout(std_out,message,'COLL')
707  endif
708 
709 
710 !=======================
711 !==  Variables for DMFT itself
712 !=======================
713 
714  mbandc = paw_dmft%mbandc
715 
716  ABI_ALLOCATE(paw_dmft%eigen_lda,(paw_dmft%nsppol,paw_dmft%nkpt,paw_dmft%mbandc))
717  paw_dmft%eigen_lda=zero
718 
719 ! allocate(paw_dmft%wtk(paw_dmft%nkpt))
720  paw_dmft%wtk=>dtset%wtk
721  if(dtset%iscf<0) then
722  paw_dmft%wtk=one/float(dtset%nkpt)
723  endif
724  sumwtk=0
725  do ikpt=1,paw_dmft%nkpt
726    sumwtk=sumwtk+paw_dmft%wtk(ikpt)
727  enddo
728  if(abs(sumwtk-1_dp)>tol13.and.dtset%iscf>=0) then
729    write(message, '(a,f15.10)' )' sum of k-point is incorrect',sumwtk
730    MSG_ERROR(message)
731  endif
732  ABI_ALLOCATE(paw_dmft%lpawu,(paw_dmft%natom))
733  do iatom=1,paw_dmft%natom
734    paw_dmft%lpawu(iatom)=pawtab(typat(iatom))%lpawu
735    if(paw_dmft%dmftqmc_t2g==1.and.paw_dmft%lpawu(iatom)==2) paw_dmft%lpawu(iatom)=1
736  enddo
737  paw_dmft%maxlpawu=maxval(paw_dmft%lpawu(:))
738  maxlpawu = paw_dmft%maxlpawu
739 
740  ABI_ALLOCATE(paw_dmft%psichi,(nsppol,dtset%nkpt,mbandc,nspinor,dtset%natom,(2*maxlpawu+1)))
741 
742  paw_dmft%psichi=cmplx(zero,zero,kind=dp)
743  paw_dmft%lpsichiortho=0
744 
745 
746 !=======================
747 ! Real      frequencies
748 !=======================
749  ABI_ALLOCATE(paw_dmft%omega_r,(2*paw_dmft%dmft_nwr))
750 
751 ! Set up real frequencies for spectral function in Hubbard one.
752  step=0.0015_dp
753  paw_dmft%omega_r(2*paw_dmft%dmft_nwr)=pi*step*(two*float(paw_dmft%dmft_nwr-1)+one)
754  do ifreq=1,2*paw_dmft%dmft_nwr-1
755   paw_dmft%omega_r(ifreq)=pi*step*(two*float(ifreq-1)+one)-paw_dmft%omega_r(2*paw_dmft%dmft_nwr)
756 !  write(std_out,*) ifreq,paw_dmft%omega_r(ifreq)
757  enddo
758 
759 
760 !=======================
761 ! Imaginary frequencies
762 !=======================
763 ! Set up log frequencies
764  if(dtset%ucrpa==0) call construct_nwlo_dmft(paw_dmft)
765 
766 
767 !=========================================================
768 !== if we use ctqmc impurity solver
769 !=========================================================
770 ! IMPORTANT : paw_dmft%hybrid is corrupted somewhere in DMFT routines on
771 ! tikal_psc and max2_open64. Use a local hybrid in qmc_prep even if not optimal.
772 ! Anyway Initializing ctqmc here is not good and produce the same result for
773 ! dmft_iter =1 which speed up the convergency ...
774 ! FIXME : Move this to init_sc_dmft and find bug
775  if(paw_dmft%dmft_solv==5) then ! CTQMC initialisation
776    write(message,'(a,2x,a,f13.5)') ch10,&
777 &  " == Initializing CTQMC"
778 !   call wrtout(std_out,message,'COLL')
779 
780    ABI_DATATYPE_ALLOCATE(paw_dmft%hybrid,(paw_dmft%natom))
781    do iatom=1,paw_dmft%natom
782      if(paw_dmft%lpawu(iatom)/=-1) then
783        nflavor=2*(2*paw_dmft%lpawu(iatom)+1)
784 #ifdef HAVE_MPI
785        call CtqmcInterface_init(paw_dmft%hybrid(iatom),paw_dmft%dmftqmc_seed,paw_dmft%dmftqmc_n, &
786 &       paw_dmft%dmftqmc_therm, paw_dmft%dmftctqmc_meas,nflavor,paw_dmft%dmftqmc_l,one/paw_dmft%temp,zero,&
787 &       std_out,paw_dmft%spacecomm)
788 #else
789        call CtqmcInterface_init(paw_dmft%hybrid(iatom),paw_dmft%dmftqmc_seed,paw_dmft%dmftqmc_n, &
790 &       paw_dmft%dmftqmc_therm, paw_dmft%dmftctqmc_meas,nflavor,paw_dmft%dmftqmc_l,one/paw_dmft%temp,zero,&
791 &       std_out)
792 #endif
793        call CtqmcInterface_setOpts(paw_dmft%hybrid(iatom),&
794                                    opt_Fk      =1,&
795 &                                  opt_order   =paw_dmft%dmftctqmc_order ,&
796 &                                  opt_movie   =paw_dmft%dmftctqmc_mov   ,&
797 &                                  opt_analysis=paw_dmft%dmftctqmc_correl,&
798 &                                  opt_check   =paw_dmft%dmftctqmc_check ,&
799 &                                  opt_noise   =paw_dmft%dmftctqmc_grnns ,&
800 &                                  opt_spectra =paw_dmft%dmftctqmc_mrka  ,&
801 &                                  opt_gmove   =paw_dmft%dmftctqmc_gmove )
802      end if
803    enddo
804    write(message,'(a,2x,a,f13.5)') ch10,&
805 &  " == Initialization CTQMC done"
806    !call wrtout(std_out,message,'COLL')
807  endif
808 
809  if(paw_dmft%dmftcheck==1.and.paw_dmft%dmft_solv<4) then
810    paw_dmft%dmftqmc_l=64
811  endif
812 
813 !************************************************************************
814 end subroutine init_dmft

m_paw_dmft/init_sc_dmft [ Functions ]

[ Top ] [ m_paw_dmft ] [ Functions ]

NAME

 init_sc_dmft

FUNCTION

  Allocate variables used in type paw_dmft_type.

INPUTS

 dmftbandi = lower bound for band states included in DMFT calculation
 dmftbandf = upper bound for band states included in DMFT calculation
 mband     = max number of bands
 nband     = number of bands for each k-point
 nkpt      = number of k-points
 nsppol    = number of spin polarisation
 occ       =  occupations
 usedmft  = ==1 if dmft is activated
 use_sc_dmft = for self-consistency in dmft

 OUTPUTS
 paw_dmft  = structure of data for dmft

PARENTS

      gstate

CHILDREN

      wrtout

SOURCE

328 subroutine init_sc_dmft(bandkss,dmftbandi,dmftbandf,dmft_read_occnd,mband,nband,nkpt,nspden,&
329 &nspinor,nsppol,occ,usedmft,paw_dmft,use_sc_dmft,dmft_solv,mpi_enreg)
330 
331 
332 !This section has been created automatically by the script Abilint (TD).
333 !Do not modify the following lines by hand.
334 #undef ABI_FUNC
335 #define ABI_FUNC 'init_sc_dmft'
336 !End of the abilint section
337 
338  implicit none
339 
340 
341 !Arguments ------------------------------------
342 !scalars
343  integer, intent(in) :: bandkss,dmft_read_occnd,dmftbandi,dmftbandf,mband,nkpt,nspden,&
344 & nspinor,nsppol,usedmft,use_sc_dmft,dmft_solv
345 !type
346  type(paw_dmft_type),intent(out) :: paw_dmft
347  type(MPI_type), intent(in) :: mpi_enreg
348 ! arrays
349  integer,intent(in) :: nband(nkpt*nsppol)
350  real(dp),intent(in) :: occ(mband*nkpt*nsppol)
351 !Local variables ------------------------------------
352  integer :: iband,icb,ikpt,isppol,nband_k,bdtot_index
353  integer :: myproc,nproc,spacecomm,use_dmft
354 ! integer :: ie,nb_procs
355  character(len=500) :: message
356 
357 !************************************************************************
358  use_dmft=abs(usedmft)
359 
360 
361 ! Check processors for DMFT
362 ! Initialise spaceComm, myproc, and nproc
363 
364  spacecomm=mpi_enreg%comm_cell
365  myproc=mpi_enreg%me_cell
366  nproc=mpi_enreg%nproc_cell
367  spacecomm=mpi_enreg%comm_world
368  myproc=mpi_enreg%me
369  nproc=mpi_enreg%nproc
370  !print *, " spacecomm,myproc,nproc",spacecomm,myproc,nproc
371  paw_dmft%spacecomm=spacecomm
372  paw_dmft%myproc=myproc
373  paw_dmft%nproc=nproc
374 
375  ! Do not comment these lines: it guarantees the parallelism in DMFT/HI or QMC will work.
376  if ((use_dmft/=0).and.(xmpi_comm_size(xmpi_world) /= xmpi_comm_size(mpi_enreg%comm_world))) then
377    MSG_ERROR("Someone changed the k-point parallelism again")
378  end if
379 
380  if(use_dmft/=0) then
381    write(message,'(2a,i3)') ch10,'-       ( number of procs used in dmft ) =', nproc
382    call wrtout(std_out,message,'COLL')
383    call wrtout(ab_out,message,'COLL')
384    write(std_out_default,'(2a,i3)') ch10,'       ( current proc is        ) =', myproc
385   ! write(ab_out_default,'(2a,i3)') ch10,'       ( current proc is        ) =', myproc
386    if(myproc==nproc-1) then
387      write(std_out_default,'(2a,i3)') ch10,'      ( last proc            ) =', myproc
388   !   write(ab_out_default,'(2a,i3)') ch10,'       ( last proc            ) =', myproc
389    endif
390  endif
391 !#ifdef HAVE_MPI
392 ! call MPI_COMM_SIZE(MPI_COMM_WORLD,nb_procs,ie)
393 ! write(6,*) "nprocs,nb_procs",nproc,nb_procs
394 ! if(nb_procs/=nproc)  then
395 !   message = ' Number of procs used in DMFT is erroneously computed '
396 !   MSG_ERROR(message)
397 ! endif
398 !#endif
399 
400  paw_dmft%mband       = mband
401  paw_dmft%dmftbandf   = dmftbandf
402  paw_dmft%dmftbandi   = dmftbandi
403  paw_dmft%nkpt        = nkpt
404 
405 !  Spin related variables and check
406  paw_dmft%nsppol      = nsppol
407  paw_dmft%nspinor     = nspinor
408  paw_dmft%nspden      = nspden
409  if(nspinor==2.and.nspden==1.and.use_dmft/=0) then
410    message = ' nspinor==2 and nspden =1 and usedmft=1 is not implemented'
411    MSG_ERROR(message)
412  endif
413 
414 ! if(nspinor==1.and.nspden==1.and.use_dmft/=0) then
415 !   message = ' nspinor==1 and nspden =1 and usedmft=1 is not implemented'
416 !   MSG_ERROR(message)
417 ! endif
418 
419  paw_dmft%use_dmft    = use_dmft
420  if (bandkss/=0) then
421    paw_dmft%use_sc_dmft = 0
422  else
423    paw_dmft%use_sc_dmft = use_sc_dmft
424  endif
425  paw_dmft%dmft_read_occnd = dmft_read_occnd
426  paw_dmft%idmftloop=0
427  paw_dmft%mbandc  = 0
428  ABI_ALLOCATE(paw_dmft%occnd,(2,mband,mband,nkpt,nsppol*use_dmft))
429  ABI_ALLOCATE(paw_dmft%band_in,(mband*use_dmft))
430  ABI_ALLOCATE(paw_dmft%include_bands,((dmftbandf-dmftbandi+1)*use_dmft))
431  ABI_ALLOCATE(paw_dmft%exclude_bands,(mband*use_dmft))
432 ! allocate(paw_dmft%ph0phiiint()
433  paw_dmft%band_in(:)=.false.
434  paw_dmft%occnd=zero
435  icb=0
436  if(use_dmft==1) then
437   do iband=1, mband
438    if(iband>=paw_dmft%dmftbandi.and.iband<=paw_dmft%dmftbandf) then
439     paw_dmft%band_in(iband)=.true.
440     paw_dmft%mbandc = paw_dmft%mbandc+1
441     paw_dmft%include_bands(paw_dmft%mbandc) = iband
442    else
443     icb=icb+1
444     paw_dmft%exclude_bands(icb)=iband
445    endif
446   enddo
447   bdtot_index=1
448   do isppol=1,nsppol
449    do ikpt=1,nkpt
450     nband_k=nband(ikpt+(isppol-1)*nkpt)
451     do iband=1,nband_k
452      paw_dmft%occnd(1,iband,iband,ikpt,isppol)=occ(bdtot_index)
453      bdtot_index=bdtot_index+1
454     end do
455    end do
456   end do
457  else
458   paw_dmft%mbandc = 0
459  endif
460  if(paw_dmft%use_dmft > 0 .and. paw_dmft%mbandc /= dmftbandf-dmftbandi+1) then
461   write(message, '(3a)' )&
462 &  ' WARNING init_sc_dmft',ch10,&
463 &  '  number of bands in dmft is not correctly computed ',ch10, &
464 &  '  Action : check the code'
465   MSG_WARNING(message)
466  endif
467  if(use_dmft>=1) then
468    write(message, '(7a)' ) ch10,ch10," ******************************************", &
469 &   ch10," DFT+DMFT Method is used", &
470 &   ch10," ******************************************"
471    call wrtout(std_out,  message,'COLL')
472    call wrtout(ab_out,  message,'COLL')
473  endif
474 
475  if(use_dmft>=1) then
476    if(dmft_solv==0) then
477      write(message, '(a,a)') ch10,' DMFT check: no solver and U=J=0'
478    else if(dmft_solv==1) then
479      write(message, '(a,a)') ch10,' DMFT check: static solver'
480    else if(dmft_solv==-1) then
481      write(message, '(a,a)') ch10,' DMFT check: static solver without renormalization of projectors: should recover LDA+U'
482    else if(dmft_solv==2) then
483      write(message, '(a,a)') ch10,' DMFT uses the Hubbard one solver'
484    else if(dmft_solv==4) then
485      write(message, '(a,a)') ch10,' DMFT uses the Hirsch Fye solver'
486    else if(dmft_solv==5) then
487      write(message, '(a,a)') ch10,' DMFT uses the Continuous Time Quantum Monte Carlo solver of ABINIT'
488    else if(dmft_solv==6) then
489      write(message, '(a,a)') ch10,' DMFT uses the Continuous Time Quantum Monte Carlo solver of TRIQS&
490      & (with density density interactions)'
491    else if(dmft_solv==7) then
492      write(message, '(a,a)') ch10,' DMFT uses the Continuous Time Quantum Monte Carlo solver of TRIQS&
493      & (with rotationaly invariant interaction)'
494   endif
495   call wrtout(std_out,message,'COLL')
496   call wrtout(ab_out,message,'COLL')
497  endif
498 
499 end subroutine init_sc_dmft

m_paw_dmft/paw_dmft_type [ Types ]

[ Top ] [ m_paw_dmft ] [ Types ]

NAME

  paw_dmft_type

FUNCTION

  This structured datatype contains the necessary data for the link
  between dmft and paw.
  occnd(non-diagonal band occupations for self-consistency), band_in
  (say which band are taken into account in the calculation), and the

SOURCE

 74  type, public :: paw_dmft_type
 75 
 76   integer :: dmft_dc
 77   ! Type of double counting used in DMFT
 78 
 79   integer :: dmft_iter
 80   ! Nb of iterations for dmft
 81 
 82   integer :: dmft_solv
 83   ! choice of solver for DMFT
 84 
 85   integer :: dmftcheck
 86   ! Check various part of the implementation
 87 
 88   integer :: dmft_log_freq
 89   ! = 0: do not use log frequencies
 90   ! = 1: use log frequencies
 91 
 92 !  integer :: dmft_mag
 93 !  ! 0 if non magnetic calculation, 1 if magnetic calculation
 94 
 95   integer :: dmft_nwlo
 96   ! dmft frequencies
 97 
 98   integer :: dmft_nwr
 99   ! dmft frequencies
100 
101   integer :: dmft_nwli
102   ! dmft frequencies
103 
104   integer :: dmftqmc_l
105   ! qmc related input
106 
107   integer :: dmftqmc_seed
108   ! qmc seed
109 
110   integer :: dmftqmc_therm
111   ! qmc thermalization
112 
113   integer :: dmftctqmc_basis
114   ! CTQMC: Basis to perform the CTQMC calculation
115   ! for historical reasons and tests
116   ! 0 : Slm basis, 1 : diagonalise local Hamiltonian, 2: diago density matrix
117 
118   integer :: dmftctqmc_check
119   ! CTQMC: perform a check on the impurity and/or bath operator
120   ! only for debug
121   ! 0 : nothing, 1 : impurity, 2 : bath, 3 : both
122 
123   integer :: dmftctqmc_correl
124   ! CTQMC: Gives analysis for CTQMC
125   ! 0 : nothing, 1 : activated Correlations.dat
126 
127   integer :: dmftctqmc_gmove
128   ! CTQMC: add global move every dmftctqmc_gmove sweeps
129   ! >= 0 ; warning inside CT-QMC with warning
130   ! == 0 ; no global moves
131 
132   integer :: dmftctqmc_grnns
133   ! CTQMC: compute green function noise for each imaginary time
134   ! 0 : nothing, 1 : activated
135 
136   integer :: dmftctqmc_meas
137   ! CTQMC: every each dmftctqmc_meas step energy is measured
138   ! Speed up caltion about 10% for dmftctqmc_meas = 2
139 
140   integer :: dmftctqmc_mrka
141   ! CTQMC: Write a temporary file Spectra_RANK.dat with the sweep evolution of
142   ! the number of electron for each flavor
143   ! The measurement is done every dmftctqmc_meas*dmftctqmc_mrka sweep
144   ! e.g. : meas=2 mrka=10 -> every 20 sweeps sum_i c+(ti)c(t'i) is mesured
145 
146   integer :: dmftctqmc_mov
147   ! CTQMC: Gives movie for CTQMC
148   ! 0 : nothing, 1 : 1 file Movie_RANK.tex for each cpu
149 
150   integer :: dmftctqmc_order
151   ! CTQMC: Gives order in perturbation for CTQMC solver
152   ! 0 : nothing, >=1 max order evaluated in Perturbation.dat
153 
154   integer :: dmftctqmc_triqs_nleg
155   ! CTQMC of TRIQS: Nb of Legendre polynomial used to compute the
156   ! Green's function (Phys. Rev. B 84, 075145) [[cite:Boehnke2011]]. Default is 30.
157   
158   ! 0 : nothing, >=1 max order evaluated in Perturbation.dat
159 
160   real(dp) :: dmftqmc_n
161   ! qmc number of sweeps
162 
163 
164   integer :: dmftqmc_t2g
165   ! for doing qmc with t2g only (for testing purposes)
166 
167   integer :: dmftbandi
168   ! Number of bands
169 
170   integer :: dmftbandf
171   ! Number of bands
172 
173   integer :: dmft_read_occnd
174   ! Number of bands
175 
176   integer :: dmft_rslf
177   ! Number of bands
178 
179   integer :: dmft_prgn
180   ! Precise the way of printing the green function.
181   !  =1   print green
182   !  =2   print self
183 
184   integer :: idmftloop
185   ! current iteration in the dmft loop
186 
187   integer :: maxlpawu         ! Number of correlated atoms
188 
189 
190   integer :: mband
191   ! Number of bands
192 
193   integer :: mbandc
194   ! Total number of bands in the Kohn-Sham Basis for PAW+DMFT
195 
196   integer :: natom
197   ! Number of atom
198 
199   integer :: natpawu         ! Number of correlated atoms
200 
201   integer :: nkpt
202   ! Number of k-point in the IBZ.
203 
204   integer :: nspden
205 
206   integer :: nspinor
207 
208   integer :: nsppol
209 
210   integer :: prtdos
211 
212   integer :: prtvol
213 
214   integer  :: lpsichiortho
215 
216   integer  :: use_fixed_self
217 
218   real(dp) :: edmft
219 
220   real(dp) :: dmft_chpr
221   ! Precision on charge required for determination of fermi level (fermi_green) with newton method
222 
223   real(dp) :: dmft_fepr
224   ! Required precision on Fermi level (fermi_green) during the DMFT SCF cycle, (=> ifermie_cv)
225   ! used also for self (new_self)  (=> iself_cv).
226 
227   real(dp) :: dmft_mxsf
228   ! Mixing coefficient for Self-Energy during the SCF DMFT cycle.
229 
230   real(dp) :: dmft_tolfreq
231   ! Required precision on local correlated density matrix  (depends on
232   ! frequency mesh), used in m_dmft/dmft_solve
233 
234   real(dp) :: dmft_lcpr
235   ! Required precision on local correlated charge  in order to stop SCF
236   ! DMFT cycle (integrate_green) => ichargeloc_cv
237 
238   real(dp) :: fermie
239 
240   real(dp) :: fermie_lda
241 
242   real(dp) :: nelectval
243 
244   character(len=fnlen) :: filapp
245 
246   real(dp) :: temp
247 
248   integer, allocatable :: lpawu(:)
249 
250   integer, allocatable :: include_bands(:)
251   ! for each bands included in the calculation (1..mbandc), include_bands
252   ! gives the index in the full band index  (1...mband)
253 
254   integer, allocatable :: exclude_bands(:)
255   ! gives the bands than are not in the DMFT calculations.
256 
257   real(dp), allocatable :: occnd(:,:,:,:,:)
258   ! non diagonal band-occupation for each k-point, polarisation.
259 
260 !  real(dp), allocatable :: phi0phiiint(:)
261 !  ! non diagonal band-occupation for each k-point, polarisation.
262 
263   logical, allocatable :: band_in(:)
264   ! true for each band included in the calculation.
265 
266   integer :: use_dmft
267   ! 1 if non diagonal occupations are used, else 0
268 
269   integer :: use_sc_dmft
270   ! 1 if calculations have to be carried out self-consistently in the
271   ! electronic density.
272 
273   complex(dpc), allocatable :: psichi(:,:,:,:,:,:)
274 
275   real(dp), allocatable :: eigen_lda(:,:,:)
276 
277   real(dp), pointer :: wtk(:) => null()
278   real(dp), pointer :: fixed_self(:,:,:,:) => null()
279   real(dp), pointer :: omega_lo(:) => null()
280   real(dp), pointer :: omega_r(:) => null()
281   real(dp), pointer :: wgt_wlo(:) => null()
282 
283   type(CtqmcInterface), allocatable :: hybrid(:)
284   type(data4entropyDMFT_t) :: forentropyDMFT
285 
286   ! MPI realated variables
287   integer :: myproc
288   integer :: nproc
289   integer :: spacecomm
290 
291  end type paw_dmft_type

m_paw_dmft/print_dmft [ Functions ]

[ Top ] [ m_paw_dmft ] [ Functions ]

NAME

 print_dmft

FUNCTION

INPUTS

OUTPUT

PARENTS

      outscfcv,vtorho

CHILDREN

      wrtout

SOURCE

1341 subroutine print_dmft(paw_dmft,pawprtvol)
1342 
1343 
1344 !This section has been created automatically by the script Abilint (TD).
1345 !Do not modify the following lines by hand.
1346 #undef ABI_FUNC
1347 #define ABI_FUNC 'print_dmft'
1348 !End of the abilint section
1349 
1350  implicit none
1351 
1352 !Arguments ------------------------------------
1353 !type
1354  type(paw_dmft_type),intent(in) :: paw_dmft
1355  integer :: pawprtvol
1356 
1357 !Local variables-------------------------------
1358  integer :: ikpt,iband,ifreq,isppol
1359  character(len=500) :: message
1360 ! *********************************************************************
1361 
1362  if( abs(pawprtvol) >= 3 )  then
1363   write(message,'(4a,3(a,2x,e21.14,a))') &
1364 &   "  -------------------------------------------------",ch10,&
1365 &   "  --- Data for DMFT ",ch10,&
1366 &   "  --- paw_dmft%fermie     = ",paw_dmft%fermie    ,ch10,&
1367 &   "  --- paw_dmft%fermie_lda = ",paw_dmft%fermie_lda,ch10,&
1368 &   "  --- paw_dmft%temp       = ",paw_dmft%temp      ,ch10
1369   call wrtout(std_out,message,'COLL')
1370   write(message,'(7(a,15x,i8,a),a,2x,e21.14,2a)') &
1371 &   "  --- paw_dmft%natpawu    = ",paw_dmft%natpawu   ,ch10,&
1372 &   "  --- paw_dmft%dmft_iter  = ",paw_dmft%dmft_iter ,ch10,&
1373 &   "  --- paw_dmft%dmft_solv  = ",paw_dmft%dmft_solv ,ch10,&
1374 &   "  --- paw_dmft%dmft_nwlo  = ",paw_dmft%dmft_nwlo ,ch10,&
1375 &   "  --- paw_dmft%dmft_nwli  = ",paw_dmft%dmft_nwli ,ch10,&
1376 &   "  --- paw_dmft%dmft_dc    = ",paw_dmft%dmft_dc   ,ch10,&
1377 &   "  --- paw_dmft%dmftqmc_l  = ",paw_dmft%dmftqmc_l ,ch10,&
1378 &   "  --- paw_dmft%dmftqmc_n  = ",paw_dmft%dmftqmc_n ,ch10,&
1379 &   "  -------------------------------------------------"
1380   call wrtout(std_out,message,'COLL')
1381 
1382 !  write(message,'(4a,3(a,2x,f8.3,a),8(a,2x,i8,a),a)') "-----------------------------------------------",ch10,&
1383 !&   "--- Data for DMFT ",ch10,&
1384 !&   "--- paw_dmft%fermie     = ",paw_dmft%fermie    ,ch10,&
1385 !&   "--- paw_dmft%fermie_lda = ",paw_dmft%fermie_lda,ch10,&
1386 !&   "--- paw_dmft%temp       = ",paw_dmft%temp      ,ch10,&
1387 !&   "--- paw_dmft%natpawu    = ",paw_dmft%natpawu   ,ch10,&
1388 !&   "--- paw_dmft%dmft_iter  = ",paw_dmft%dmft_iter ,ch10,&
1389 !&   "--- paw_dmft%dmft_solv  = ",paw_dmft%dmft_solv ,ch10,&
1390 !&   "--- paw_dmft%dmft_nwlo  = ",paw_dmft%dmft_nwlo ,ch10,&
1391 !&   "--- paw_dmft%dmft_nwli  = ",paw_dmft%dmft_nwli ,ch10,&
1392 !&   "--- paw_dmft%dmft_dc    = ",paw_dmft%dmft_dc   ,ch10,&
1393 !&   "--- paw_dmft%dmftqmc_l  = ",paw_dmft%dmftqmc_l ,ch10,&
1394 !&   "--- paw_dmft%dmftqmc_n  = ",paw_dmft%dmftqmc_n ,ch10,&
1395 !&   "-----------------------------------------------"
1396   if(abs(pawprtvol)>10) then
1397    call wrtout(std_out,message,'COLL')
1398    write(message, '(a)') " LDA Eigenvalues "
1399    do isppol=1,paw_dmft%nsppol
1400     write(message, '(a,i4)') "--isppol--",isppol
1401     call wrtout(std_out,message,'COLL')
1402     do ikpt=1,paw_dmft%nkpt
1403      write(message, '(a,i4,2x,f14.5,a)') "  -k-pt--",ikpt,paw_dmft%wtk(ikpt),"(<-weight(k-pt))"
1404 
1405      call wrtout(std_out,message,'COLL')
1406      do iband=1,paw_dmft%mbandc
1407       write(message, '(a,i4,f10.5)') "   -iband--",iband,paw_dmft%eigen_lda(isppol,ikpt,iband)
1408       call wrtout(std_out,message,'COLL')
1409      enddo
1410     enddo
1411    enddo
1412    write(message, '(3x,a)') "Log. Freq"
1413    call wrtout(std_out,message,'COLL')
1414    do ifreq=1,paw_dmft%dmft_nwlo
1415     write(message, '(3x,a,i4,2(2x,e13.5))') "--ifreq--",ifreq,paw_dmft%omega_lo(ifreq),paw_dmft%wgt_wlo(ifreq)
1416     call wrtout(std_out,message,'COLL')
1417    enddo
1418   endif
1419  endif
1420 
1421 end subroutine print_dmft

m_paw_dmft/print_sc_dmft [ Functions ]

[ Top ] [ m_paw_dmft ] [ Functions ]

NAME

 print_sc_dmft

FUNCTION

INPUTS

OUTPUT

PARENTS

      gstate

CHILDREN

      wrtout

SOURCE

1442 subroutine print_sc_dmft(paw_dmft,pawprtvol)
1443 
1444 
1445 !This section has been created automatically by the script Abilint (TD).
1446 !Do not modify the following lines by hand.
1447 #undef ABI_FUNC
1448 #define ABI_FUNC 'print_sc_dmft'
1449 !End of the abilint section
1450 
1451  implicit none
1452 
1453 !Arguments ------------------------------------
1454 !type
1455  type(paw_dmft_type),intent(in) :: paw_dmft
1456  integer :: pawprtvol
1457 
1458 !Local variables-------------------------------
1459  integer :: iband
1460  character(len=500) :: message
1461 ! *********************************************************************
1462 
1463  if( abs(pawprtvol) >= 3 )  then
1464    write(message,'(5a,8(a,2x,i5,a),a)')ch10,"-----------------------------------------------",ch10,&
1465 &    "--- Data for SC DMFT ",ch10,&
1466 &    "--- paw_dmft%mband       = ",paw_dmft%mband,ch10,&
1467 &    "--- paw_dmft%dmftbandf   = ",paw_dmft%dmftbandf,ch10,&
1468 &    "--- paw_dmft%dmftbandi   = ",paw_dmft%dmftbandi,ch10,&
1469 &    "--- paw_dmft%nkpt        = ",paw_dmft%nkpt,ch10,&
1470 &    "--- paw_dmft%nsppol      = ",paw_dmft%nsppol,ch10,&
1471 &    "--- paw_dmft%use_dmft    = ",paw_dmft%use_dmft,ch10,&
1472 &    "--- paw_dmft%use_sc_dmft = ",paw_dmft%use_sc_dmft,ch10,&
1473 &    "--- paw_dmft%mbandc      = ",paw_dmft%mbandc,ch10,&
1474 &    "-----------------------------------------------"
1475    call wrtout(std_out,message,'COLL')
1476    write(message, '(a)') " paw_dmft%band_in"
1477    call wrtout(std_out,message,'COLL')
1478    write(message, '(100i5)') (iband,iband=1,min(paw_dmft%mband,100))
1479    call wrtout(std_out,message,'COLL')
1480    write(message, '(100L3)') (paw_dmft%band_in(iband),iband=1,min(paw_dmft%mband,100))
1481    call wrtout(std_out,message,'COLL')
1482    do iband=1,paw_dmft%mbandc
1483      write(message,*) "include_bands",iband,paw_dmft%include_bands(iband)
1484      call wrtout(std_out,message,'COLL')
1485    enddo
1486    write(message, '(a,a,i4,a)' )ch10,&
1487 &    'The',paw_dmft%mband-paw_dmft%dmftbandf+paw_dmft%dmftbandi-1,&
1488 &    '  Following bands are excluded from the DMFT calculation  '
1489    call wrtout(std_out,message,'COLL')
1490    write(message,'(100i5)') (paw_dmft%exclude_bands(iband),iband=1,min(paw_dmft%mband-paw_dmft%dmftbandf+paw_dmft%dmftbandi-1,100))
1491    call wrtout(std_out,message,'COLL')
1492  endif
1493 
1494 end subroutine print_sc_dmft

m_paw_dmft/readocc_dmft [ Functions ]

[ Top ] [ m_paw_dmft ] [ Functions ]

NAME

 readocc_dmft

FUNCTION

  read occnd on disk

INPUTS

  paw_dmft   = data structure
  filnam_ds3 = root for filname to read (input)
  filnam_ds4 = root for filname to read (output)

OUTPUT

  paw_dmft: occnd

PARENTS

      gstate

CHILDREN

      wrtout

SOURCE

1594 subroutine readocc_dmft(paw_dmft,filnam_ds3,filnam_ds4)
1595 
1596 
1597 !This section has been created automatically by the script Abilint (TD).
1598 !Do not modify the following lines by hand.
1599 #undef ABI_FUNC
1600 #define ABI_FUNC 'readocc_dmft'
1601 !End of the abilint section
1602 
1603  implicit none
1604 
1605 !Arguments ------------------------------------
1606  type(paw_dmft_type),intent(inout) :: paw_dmft
1607  character(len=fnlen) :: filnam_ds3
1608  character(len=fnlen) :: filnam_ds4
1609 
1610 !Local variables-------------------------------
1611 !scalars
1612  character(len=fnlen) :: tmpfil
1613  integer :: ib,ib1,ikpt,is,unitsaveocc,dum1,dum2,dum3,dum4,ioerr
1614  logical :: lexist
1615  character(len=500) :: message
1616  character(len=4) :: chtemp
1617 ! *********************************************************************
1618  if(paw_dmft%dmft_read_occnd==0) return
1619  if(paw_dmft%dmft_read_occnd==1) tmpfil=trim(filnam_ds3)//'_DMFTOCCND'
1620  if(paw_dmft%dmft_read_occnd==2) tmpfil=trim(filnam_ds4)//'_DMFTOCCND'
1621  inquire(file=trim(tmpfil),exist=lexist)!,recl=nrecl)
1622  unitsaveocc=679
1623  if (lexist) then
1624    if (open_file(tmpfil,message,unit=unitsaveocc,status='unknown',form='formatted') /= 0) then
1625      MSG_ERROR(message)
1626    end if
1627    rewind(unitsaveocc)
1628    write(message,'(3a)') ch10,"  == Read DMFT non diagonal occupations on disk"
1629    call wrtout(std_out,message,'COLL')
1630    read(unitsaveocc,*)
1631    read(unitsaveocc,*,iostat=ioerr)&
1632 &              chtemp,dum1,dum2,dum3,dum4
1633    if(ioerr<0) then
1634      write(std_out,*) "read",dum1,dum2,dum3,dum4
1635    endif
1636    write(message,'(2a,4i4)') ch10,"  == natom, nsppol, nbandc, nkpt  read are",dum1,dum2,dum3,dum4
1637    call wrtout(std_out,message,'COLL')
1638    do is = 1 , paw_dmft%nsppol
1639      do ikpt = 1, paw_dmft%nkpt
1640        do ib = 1, paw_dmft%mbandc
1641          do ib1 = 1, paw_dmft%mbandc
1642            if (paw_dmft%nspinor==1) then
1643              read(unitsaveocc,*) dum1,dum2,dum3,dum4,&
1644 &             paw_dmft%occnd(1,paw_dmft%include_bands(ib),paw_dmft%include_bands(ib1),ikpt,is)
1645            endif
1646            if (paw_dmft%nspinor==2) then
1647              read(unitsaveocc,*) dum1,dum2,dum3,dum4,&
1648 &             paw_dmft%occnd(1,paw_dmft%include_bands(ib),paw_dmft%include_bands(ib1),ikpt,is),&
1649 &             paw_dmft%occnd(2,paw_dmft%include_bands(ib),paw_dmft%include_bands(ib1),ikpt,is)
1650            endif
1651          enddo
1652        enddo
1653      enddo
1654    enddo
1655 !   write(read,'(3a)') "# end of record",ch10&
1656 !&                ,"####  1234 "
1657 !   call wrtout(unitsaveocc,message,'COLL')
1658  else
1659    write(message,'(2a,2x,2a)') ch10,"   File",trim(tmpfil),"is not available"
1660    call wrtout(std_out,message,'COLL')
1661    write(message,'(4a)') ch10,"  ==> DMFT Occupations not available for restart", ch10, &
1662 &   "      -> The calculation is started with Fermi Dirac scheme for occupations"
1663    call wrtout(std_out,message,'COLL')
1664  endif
1665 
1666 end subroutine readocc_dmft

m_paw_dmft/saveocc_dmft [ Functions ]

[ Top ] [ m_paw_dmft ] [ Functions ]

NAME

 saveocc_dmft

FUNCTION

  save occnd on disk

INPUTS

  paw_dmft

OUTPUT

PARENTS

      vtorho

CHILDREN

      wrtout

SOURCE

1517 subroutine saveocc_dmft(paw_dmft)
1518 
1519 
1520 !This section has been created automatically by the script Abilint (TD).
1521 !Do not modify the following lines by hand.
1522 #undef ABI_FUNC
1523 #define ABI_FUNC 'saveocc_dmft'
1524 !End of the abilint section
1525 
1526  implicit none
1527 
1528 !Arguments ------------------------------------
1529  type(paw_dmft_type),intent(inout) :: paw_dmft
1530 
1531 !Local variables-------------------------------
1532 !scalars
1533  character(len=fnlen) :: tmpfil
1534  integer :: ib,ib1,ikpt,is,unitsaveocc
1535  character(len=500) :: message
1536 ! *********************************************************************
1537  tmpfil = trim(paw_dmft%filapp)//'_DMFTOCCND'
1538  if (open_file(tmpfil,message,newunit=unitsaveocc,status='unknown',form='formatted') /= 0) then
1539    MSG_ERROR(message)
1540  end if
1541 
1542  rewind(unitsaveocc)
1543  write(message,'(2a)') ch10,"  == Print DMFT non diagonal occupations on disk"
1544  call wrtout(std_out,message,'COLL')
1545  write(message,'(3a,2x,4i5)') "# natom,nsppol,mbandc,nkpt",ch10&
1546 &              ,"####",paw_dmft%natom,paw_dmft%nsppol,paw_dmft%mbandc,paw_dmft%nkpt
1547  call wrtout(unitsaveocc,message,'COLL')
1548  do is = 1 , paw_dmft%nsppol
1549    do ikpt = 1, paw_dmft%nkpt
1550      do ib = 1, paw_dmft%mbandc
1551        do ib1 = 1, paw_dmft%mbandc
1552          if (paw_dmft%nspinor==1) then
1553            write(unitsaveocc,*) is,ikpt,ib,ib1,paw_dmft%occnd(1,paw_dmft%include_bands(ib),paw_dmft%include_bands(ib1),ikpt,is)
1554          endif
1555          if (paw_dmft%nspinor==2) then
1556            write(unitsaveocc,*) is,ikpt,ib,ib1,paw_dmft%occnd(1,paw_dmft%include_bands(ib),paw_dmft%include_bands(ib1),ikpt,is),&
1557 &           paw_dmft%occnd(2,paw_dmft%include_bands(ib),paw_dmft%include_bands(ib1),ikpt,is)
1558          endif
1559        enddo
1560      enddo
1561    enddo
1562  enddo
1563  write(message,'(3a)') "# end of record",ch10&
1564 &              ,"####  1234 "
1565  call wrtout(unitsaveocc,message,'COLL')
1566  close(unitsaveocc)
1567 
1568 end subroutine saveocc_dmft