TABLE OF CONTENTS


ABINIT/m_entropyDMFT [ Modules ]

[ Top ] [ Modules ]

NAME

  m_entropyDMFT

FUNCTION

  FIXME: add description.

COPYRIGHT

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

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

25 #if defined HAVE_CONFIG_H
26 #include "config.h"
27 #endif
28 
29 #include "abi_common.h"
30 
31 module m_entropyDMFT
32 
33   use defs_basis
34   use defs_abitypes
35   use m_errors
36   use m_xmpi
37   use m_energies, only : energies_type, energies_eval_eint
38   use m_splines, only : spline_integrate, spline, splint
39   use m_pawang, only : pawang_type
40   use m_pawrad, only : pawrad_type, simp_gen, poisson
41   use m_pawtab, only : pawtab_type
42   use m_paw_correlations,only : pawpuxinit
43   use m_io_tools, only : get_unit
44   use m_data4entropyDMFT
45 
46   implicit none
47 
48   private
49 
50   public :: entropyDMFT_init
51   public :: entropyDMFT_destroy
52   public :: entropyDMFT_nextLambda
53   public :: entropyDMFT_addIntegrand
54   public :: entropyDMFT_computeEntropy
55 
56   integer, parameter :: E_U0       = 1
57   integer, parameter :: E_UU       = 2
58   integer, parameter :: E_DIRECT   = 1
59   integer, parameter :: E_DC       = 2
60   integer, parameter :: AC_NOTHING = 0
61   integer, parameter :: AC_ETOT    = 1
62   character(len=21), parameter :: HDR_NAME = "DATA FOR ETOT DMFT v="

ABINIT/m_entropyDMFT/entropyDMFT_addIntegrand [ Functions ]

[ Top ] [ Functions ]

NAME

  entropyDMFT_addIntegrand

FUNCTION

  FIXME: add description.

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (J. Bieder)
  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

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

      m_scfcv

CHILDREN

SOURCE

793   subroutine entropyDMFT_addIntegrand(e_t,dt,energies,data4etot)
794 
795 !Arguments ------------------------------------
796 
797 !This section has been created automatically by the script Abilint (TD).
798 !Do not modify the following lines by hand.
799 #undef ABI_FUNC
800 #define ABI_FUNC 'entropyDMFT_addIntegrand'
801 !End of the abilint section
802 
803     type(entropyDMFT_t)   , intent(inout) :: e_t
804     type(dataset_type) , intent(in   ) :: dt
805     type(energies_type), intent(in   ) :: energies
806     type(data4entropyDMFT_t)  , intent(in   ) :: data4etot
807 !Local variables ------------------------------
808     integer :: optdc
809     integer :: ictypat
810     integer :: iatom, icatom
811     integer :: iflavor1, iflavor2, ndim
812     integer :: icouple
813 
814     if ( e_t%action == AC_NOTHING ) return
815 
816     ! Write lambda for restart
817     if ( e_t%rank == 0 ) then
818       open(unit=e_t%ofile,file=e_t%filedata,position="append",form="unformatted")
819       write(e_t%ofile) e_t%lambda(e_t%mylambda)
820       close(e_t%ofile)
821     end if
822 
823     if ( e_t%mylambda == 1 ) then
824       ! Save entropy and internal energy for U=0
825       e_t%entropy0 = energies%entropy
826       ! 1 is for usepaw that is 1 in DMFT, optdc is to know if the DC scheme is
827       ! calculated.
828       call energies_eval_eint(energies,dt,1,optdc,e_t%energies(E_DIRECT,E_U0),e_t%energies(E_DC,E_U0))
829       if ( e_t%rank == 0 ) then
830         open(unit=e_t%ofile,file=e_t%filename,position="append")
831         write(e_t%ofile,'(a)') "# Temperature [Ha]:"
832         write(e_t%ofile,'(es22.14)') e_t%temp
833         write(e_t%ofile,'(a)') "# Entropy for lambda=0 [kb]:"
834         write(e_t%ofile,'(es22.14)') e_t%entropy0
835         write(e_t%ofile,'(a)') "# Internal energy for lambda=0 [Ha]:"
836         write(e_t%ofile,'(es22.14)') e_t%energies(E_DC,E_U0)
837         close(e_t%ofile)
838         open(unit=e_t%ofile,file=e_t%filedata,position="append",form="unformatted")
839         write(e_t%ofile) e_t%temp
840         write(e_t%ofile) e_t%entropy0
841         write(e_t%ofile) e_t%energies(E_DC,E_U0)
842         close(e_t%ofile)
843       end if
844     else if ( e_t%mylambda == e_t%nlambda ) then
845       ! Save internal energy for U=Umax
846       call energies_eval_eint(energies,dt,1,optdc,e_t%energies(E_DIRECT,E_UU),e_t%energies(E_DC,E_UU))
847       if ( e_t%rank == 0 ) then
848         open(unit=e_t%ofile,file=e_t%filename,position="append")
849         write(e_t%ofile,'(a)') "# Internal energy for lambda=1 [Ha]:"
850         write(e_t%ofile,'(es22.14)') e_t%energies(E_DC,E_UU)
851       end if
852       ! Compute U_max,  and uij
853       do ictypat = 1, e_t%nctypat
854         ndim = 2*(2*e_t%lpawu(ictypat)+1)
855         icouple = 0
856         !write(*,*) matU
857         if ( e_t%rank == 0 ) &
858            write(e_t%ofile,'(a,f7.5,1x,a,i4)') "# Interaction Matrix normalized by U=",e_t%U_input(ictypat) , &
859            "[Ha] for atom type", e_t%index_typat(ictypat)
860         do iflavor1 = 1, ndim
861           if ( e_t%rank == 0 ) &
862             write(e_t%ofile,'(14(a21,2x))',advance="no") (/ ( "...", iflavor2=1,iflavor1 ) /)
863           do iflavor2 = iflavor1+1, ndim
864             icouple = icouple + 1
865             e_t%uij(icouple,ictypat) = data4etot%hu_dens(iflavor1,iflavor2,e_t%index_typat(ictypat)) &
866                                       /e_t%U_input(ictypat)   ! to modify, this is just the idea
867             if ( e_t%rank == 0 ) &
868               write(e_t%ofile,'(14(es21.14,2x))',advance="no") e_t%uij(icouple,ictypat)
869           end do
870           if ( e_t%rank == 0 ) write(e_t%ofile,*)
871         end do
872       end do
873 
874       if ( e_t%rank == 0 ) then
875         close(e_t%ofile)
876         open(unit=e_t%ofile,file=e_t%filedata,position="append",form="unformatted")
877         write(e_t%ofile) e_t%energies(E_DC,E_UU)
878         do ictypat = 1, e_t%nctypat
879           ndim = 2*(2*e_t%lpawu(ictypat)+1)
880           icouple = 0
881           do iflavor1 = 1, ndim
882             do iflavor2 = iflavor1+1, ndim
883               icouple = icouple + 1
884               write(e_t%ofile) e_t%uij(icouple,ictypat)
885             end do
886           end do
887         end do
888         close(e_t%ofile)
889       end if
890     endif
891 
892     ! For all lambda
893     ! Save Docc, Nup and Ndwn
894     if ( e_t%rank == 0 ) then
895       open(unit=e_t%ofile,file=e_t%filedata,position="append",form="unformatted")
896     end if
897     do icatom = 1, e_t%ncatom
898       iatom = e_t%index_atom(icatom)
899       ndim = 2*(2*e_t%lpawu(e_t%typat(icatom))+1)
900       icouple = 0
901       e_t%e_dc(icatom,e_t%mylambda) = data4etot%e_dc(iatom)
902       if ( e_t%rank == 0 ) then
903         write(e_t%ofile) e_t%e_dc(icatom,e_t%mylambda)
904       end if
905       do iflavor1 = 1, ndim
906         do iflavor2 = iflavor1+1, ndim
907           icouple = icouple + 1
908           e_t%docc(icouple,icatom,e_t%mylambda) = data4etot%Docc(iflavor1,iflavor2,iatom)
909           if ( e_t%rank == 0 ) then
910             write(e_t%ofile) e_t%docc(icouple,icatom,e_t%mylambda)
911           end if
912         end do
913       end do
914     end do
915 
916     if ( e_t%rank == 0 ) then
917       close(e_t%ofile)
918     end if
919 
920     !call entropyDMFT_computeIntegrand(e_t)
921 
922   end subroutine entropyDMFT_addIntegrand

ABINIT/m_entropyDMFT/entropyDMFT_allocateAll [ Functions ]

[ Top ] [ Functions ]

NAME

  entropyDMFT_allocateAll

FUNCTION

  FIXME: add description.

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (J. Bieder)
  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

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

      m_entropyDMFT

CHILDREN

SOURCE

339   subroutine entropyDMFT_allocateAll(e_t)
340 
341 !Arguments ------------------------------------
342 
343 !This section has been created automatically by the script Abilint (TD).
344 !Do not modify the following lines by hand.
345 #undef ABI_FUNC
346 #define ABI_FUNC 'entropyDMFT_allocateAll'
347 !End of the abilint section
348 
349     type(entropyDMFT_t), intent(inout) :: e_t
350 
351     ABI_ALLOCATE(e_t%index_atom, (1:e_t%natom))
352     e_t%index_atom = 0
353     ABI_ALLOCATE(e_t%index_typat,(1:e_t%nctypat))
354     e_t%index_atom = 0
355     ABI_ALLOCATE(e_t%typat,      (1:e_t%ncatom))
356     e_t%typat      = 0
357     ABI_ALLOCATE(e_t%lpawu,      (1:e_t%nctypat))
358     e_t%lpawu      = 0
359     ABI_ALLOCATE(e_t%U_input,    (1:e_t%nctypat))
360     e_t%U_input    = zero
361     ABI_ALLOCATE(e_t%J_input,    (1:e_t%nctypat))
362     e_t%J_input    = zero
363     ABI_ALLOCATE(e_t%lambda,     (1:e_t%nlambda))
364     e_t%lambda     = 0
365     ABI_ALLOCATE(e_t%docc,       (1:(14*13)/2,1:e_t%ncatom,1:e_t%nlambda))
366     e_t%docc       = zero
367     ABI_ALLOCATE(e_t%e_dc,       (1:e_t%ncatom,1:e_t%nlambda))
368     e_t%e_dc       = zero
369     ABI_ALLOCATE(e_t%uij,        (1:(14*13)/2,1:e_t%nctypat))
370     e_t%uij        = zero
371   end subroutine entropyDMFT_allocateAll

ABINIT/m_entropyDMFT/entropyDMFT_computeEntropy [ Functions ]

[ Top ] [ Functions ]

NAME

  entropyDMFT_computeEntropy

FUNCTION

  FIXME: add description.

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (J. Bieder)
  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

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

      m_scfcv

CHILDREN

SOURCE

 955   subroutine entropyDMFT_computeEntropy(e_t,entropy)
 956 
 957 !Arguments ------------------------------------
 958 
 959 !This section has been created automatically by the script Abilint (TD).
 960 !Do not modify the following lines by hand.
 961 #undef ABI_FUNC
 962 #define ABI_FUNC 'entropyDMFT_computeEntropy'
 963 !End of the abilint section
 964 
 965     type(entropyDMFT_t), intent(inout) :: e_t
 966     real(dp)        , intent(inout) :: entropy !inout in case we do nothing, avoid to change entropy in NaN
 967 !Local variables ------------------------------
 968     integer :: ilambda
 969     integer :: nlambda
 970     integer :: icatom
 971     integer :: ncatom
 972     integer :: icouple
 973     integer :: ndim
 974     real(dp), allocatable :: integrand(:,:)
 975     real(dp) :: docc
 976     real(dp) :: integral
 977     real(dp) :: entropyDirect
 978     character(len=500) :: string
 979     character(len=50) :: i2str
 980 
 981     if ( e_t%action == AC_NOTHING ) return
 982 
 983     nlambda = e_t%nlambda
 984     ncatom  = e_t%ncatom
 985 
 986     ABI_ALLOCATE(integrand,(1:nlambda,1:ncatom))
 987     integrand(1:nlambda,1:ncatom) = zero
 988 
 989     if ( e_t%rank == 0 ) then
 990       open(unit=e_t%ofile,file=e_t%filename,position="append")
 991       write(e_t%ofile,'(2a)') ch10,"#Decomposition for each lambda per atom"
 992     end if
 993     do ilambda = 1, nlambda
 994       if ( e_t%rank == 0 ) then
 995         write(e_t%ofile,'(a,i2,a)') "# ------| lambda ",ilambda," |------ #"
 996         write(e_t%ofile,'(a5,a23,a23)') "#Atom","DC","<uij*n_i*n_j>"
 997       end if
 998       do icatom = 1, ncatom
 999         integrand(ilambda,icatom) = -e_t%e_dc(icatom,ilambda)
1000         ndim = 2*(2*e_t%lpawu(e_t%typat(e_t%index_atom(icatom)))+1) ! nflavors
1001         ndim = ndim*(ndim-1)/2 ! number of couples
1002         docc = zero
1003         do icouple = 1, ndim
1004           docc = docc + e_t%uij(icouple,e_t%typat(icatom)) * e_t%docc(icouple,icatom,ilambda)
1005         end do
1006         integrand(ilambda,icatom) = integrand(ilambda,icatom) + docc
1007         if ( e_t%rank == 0 ) &
1008           write(e_t%ofile,'(i5,2es23.14)') e_t%index_atom(icatom),e_t%e_dc(icatom,ilambda),docc
1009       end do
1010     end do
1011 
1012     ! print for test purpose
1013     if ( e_t%rank == 0 ) then
1014       write(i2str,'(i6)') e_t%ncatom
1015       string='(a,a7,'//TRIM(ADJUSTL(i2str))//'(12x,"Atom",i6))'
1016       !open(unit=e_t%ofile,file=e_t%filename,position="append")
1017       write(e_t%ofile,string) ch10,"#lambda", (/ (e_t%index_atom(icatom),icatom=1,ncatom) /)
1018       string='(1x,f6.4,'//TRIM(ADJUSTL(i2str))//'(1x,e21.14))'
1019       do ilambda = 1, nlambda
1020         write(e_t%ofile,string) e_t%lambda(ilambda), (/ (integrand(ilambda,icatom),icatom=1,ncatom) /)
1021       end do
1022       close(e_t%ofile)
1023     end if
1024 
1025     ! Integrate the integrand for all correlated atoms
1026     call entropyDMFT_integrate(e_t,integrand,integral)
1027     ABI_DEALLOCATE(integrand)
1028 
1029     write(string,'(a,1x,78a)') ch10,"+",(/ ("-",ilambda=1,76) /), "+"
1030     call wrtout(std_out,string,"COLL")
1031     call wrtout(ab_out,string,"COLL")
1032     write(string,'(1x,a)') "|             Calculation of entropy within  the DMFT Framework              |"
1033     call wrtout(std_out,string,"COLL")
1034     call wrtout(ab_out,string,"COLL")
1035     write(string,'(1x,40a)') "+",(/ ("- ",ilambda=1,38) /), "+"
1036     call wrtout(std_out,string,"COLL")
1037 
1038     entropyDirect = e_t%entropy0 + ( e_t%energies(E_DIRECT,E_UU) - e_t%energies(E_DIRECT,E_U0) - integral ) /e_t%temp
1039     entropy = e_t%entropy0 + ( e_t%energies(E_DC,E_UU) - e_t%energies(E_DC,E_U0) - integral ) /e_t%temp
1040 
1041     write(i2str,'(a)') '(1x,a,19x,a11,es21.14,1x,a4,20x,a)'
1042     write(string,i2str) "|","Integral = ", integral, "[Ha]", "|"
1043     call wrtout(std_out,string,"COLL")
1044     write(string,i2str) "|","E(0)     = ", e_t%energies(E_DC,E_U0), "[Ha]", "|"
1045     call wrtout(std_out,string,"COLL")
1046     write(string,i2str) "|","E(U)     = ", e_t%energies(E_DC,E_UU), "[Ha]", "|"
1047     call wrtout(std_out,string,"COLL")
1048     write(string,i2str) "|","S(0)     = ", e_t%entropy0, "[kb]", "|"
1049     call wrtout(std_out,string,"COLL")
1050     write(string,i2str) "|","S(U)     = ", entropy, "[kb]", "|"
1051     call wrtout(std_out,string,"COLL")
1052     write(string,'(1x,40a)') "+",(/ ("- ",ilambda=1,38) /), "+"
1053     call wrtout(std_out,string,"COLL")
1054 
1055     write(string,'(1x,a,16x,a22,es21.14,1x,a4,12x,a)') "|","-kT*entropy is set to ", -e_t%temp*entropy, "[Ha]", "|"
1056     call wrtout(std_out,string,"COLL")
1057     call wrtout(ab_out,string,"COLL")
1058 
1059     write(string,'(1x,78a)') "+",(/ ("-",ilambda=1,76) /), "+"
1060     call wrtout(std_out,string,"COLL")
1061     call wrtout(ab_out,string,"COLL")
1062 
1063     if ( entropy < zero ) then
1064       write(string,'(3a)') "Entropy is negative !!!!",ch10,&
1065       "It does not make any sense"
1066       MSG_WARNING(string)
1067     end if
1068 
1069     if ( abs(entropy-entropyDirect) >= tol3 ) then
1070       write(string,'(1x,a,1x,f8.3,1x,a,1x,f8.3,2a)') "Difference between Direct and DC entropies is", abs(entropy-entropyDirect), &
1071         "which is greater than", tol3,ch10,"Action : converge better the DMFT and/or DFT loops"
1072       MSG_WARNING(string)
1073     end if
1074 
1075 
1076 
1077   end subroutine entropyDMFT_computeEntropy

ABINIT/m_entropyDMFT/entropyDMFT_destroy [ Functions ]

[ Top ] [ Functions ]

NAME

  entropyDMFT_destroy

FUNCTION

  FIXME: add description.

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (J. Bieder)
  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

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

      m_scfcv

CHILDREN

SOURCE

1213   subroutine entropyDMFT_destroy(e_t)
1214 
1215 !Arguments ------------------------------------
1216 
1217 !This section has been created automatically by the script Abilint (TD).
1218 !Do not modify the following lines by hand.
1219 #undef ABI_FUNC
1220 #define ABI_FUNC 'entropyDMFT_destroy'
1221 !End of the abilint section
1222 
1223     type(entropyDMFT_t), intent(inout) :: e_t
1224 
1225     if ( allocated(e_t%index_atom) ) then
1226       ABI_DEALLOCATE(e_t%index_atom)
1227     endif
1228     if ( allocated(e_t%index_typat)) then
1229       ABI_DEALLOCATE(e_t%index_typat)
1230     endif
1231     if ( allocated(e_t%typat     ) ) then
1232       ABI_DEALLOCATE(e_t%typat     )
1233     endif
1234     if ( allocated(e_t%lpawu     ) ) then
1235       ABI_DEALLOCATE(e_t%lpawu)
1236     endif
1237     if ( allocated(e_t%U_input   ) ) then
1238       ABI_DEALLOCATE(e_t%U_input)
1239     endif
1240     if ( allocated(e_t%J_input   ) ) then
1241       ABI_DEALLOCATE(e_t%J_input)
1242     endif
1243     if ( allocated(e_t%lambda    ) ) then
1244       ABI_DEALLOCATE(e_t%lambda)
1245     endif
1246     if ( allocated(e_t%docc      ) ) then
1247       ABI_DEALLOCATE(e_t%docc)
1248     endif
1249     if ( allocated(e_t%e_dc      ) ) then
1250       ABI_DEALLOCATE(e_t%e_dc)
1251     endif
1252     if ( allocated(e_t%uij       ) ) then
1253       ABI_DEALLOCATE(e_t%uij)
1254     endif
1255     e_t%isset = .FALSE.
1256   end subroutine entropyDMFT_destroy

ABINIT/m_entropyDMFT/entropyDMFT_dump [ Functions ]

[ Top ] [ Functions ]

NAME

  entropyDMFT_dump

FUNCTION

  FIXME: add description.

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (J. Bieder)
  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

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

      m_entropyDMFT

CHILDREN

SOURCE

558   subroutine entropyDMFT_dump(e_t)
559 
560 !Arguments ------------------------------------
561 
562 !This section has been created automatically by the script Abilint (TD).
563 !Do not modify the following lines by hand.
564 #undef ABI_FUNC
565 #define ABI_FUNC 'entropyDMFT_dump'
566 !End of the abilint section
567 
568     type(entropyDMFT_t), intent(inout) :: e_t
569 !Local variables-------------------------------
570     integer            :: ilambda
571     integer            :: ictypat
572     integer            :: ndim
573     integer            :: icouple
574     integer            :: iflavor1
575     integer            :: iflavor2
576     integer            :: icatom
577     integer            :: iatom
578 
579     if ( e_t%rank /= 0 ) return
580 
581     ! Dump _data4EtotDMFT
582     open(unit=e_t%ofile,file=e_t%filedata,action="write",form="unformatted")
583     write(e_t%ofile) HDR_NAME, 1
584 
585     do ilambda = 1, e_t%mylambda
586       write(e_t%ofile) e_t%lambda(ilambda)
587 
588       if ( ilambda == 1 ) then
589         write(e_t%ofile) e_t%temp
590         write(e_t%ofile) e_t%entropy0
591         write(e_t%ofile) e_t%energies(E_DC,E_U0)
592       else if ( ilambda == e_t%nlambda ) then ! should never happend ?!
593         write(e_t%ofile) e_t%energies(E_DC,E_UU)
594         do ictypat = 1, e_t%nctypat
595           ndim = 2*(2*e_t%lpawu(ictypat)+1)
596           icouple = 0
597           do iflavor1 = 1, ndim
598             do iflavor2 = iflavor1+1, ndim
599               icouple = icouple + 1
600               write(e_t%ofile) e_t%uij(icouple,ictypat)
601             end do
602           end do
603         end do
604       end if
605 
606       do icatom = 1, e_t%ncatom
607         iatom = e_t%index_atom(icatom)
608         ndim = 2*(2*e_t%lpawu(e_t%typat(icatom))+1)
609         icouple = 0
610         write(e_t%ofile) e_t%e_dc(icatom,ilambda)
611         do iflavor1 = 1, ndim
612           do iflavor2 = iflavor1+1, ndim
613             icouple = icouple + 1
614             write(e_t%ofile) e_t%docc(icouple,icatom,ilambda)
615           end do
616         end do
617       end do
618     end do
619     close(e_t%ofile)
620 
621     ! Dump _EtotDMFT
622     open(unit=e_t%ofile,file=e_t%filename)
623     write(e_t%ofile,'(2a)') "# Data for entropy calculation in DMFT",ch10
624     do ilambda = 1, e_t%mylambda
625       if ( ilambda == 1 ) then
626         write(e_t%ofile,'(a)') "# Temperature [Ha]:"
627         write(e_t%ofile,'(es22.14)') e_t%temp
628         write(e_t%ofile,'(a)') "# Entropy for lambda=0 [kb]:"
629         write(e_t%ofile,'(es22.14)') e_t%entropy0
630         write(e_t%ofile,'(a)') "# Internal energy for lambda=0 [Ha]:"
631         write(e_t%ofile,'(es22.14)') e_t%energies(E_DC,E_U0)
632       else if ( e_t%mylambda == e_t%nlambda ) then
633         write(e_t%ofile,'(a)') "# Internal energy for lambda=1 [Ha]:"
634         write(e_t%ofile,'(es22.14)') e_t%energies(E_DC,E_UU)
635         do ictypat = 1, e_t%nctypat
636           ndim = 2*(2*e_t%lpawu(ictypat)+1)
637           icouple = 0
638           write(e_t%ofile,'(a,f7.5,1x,a,i4)') "# Interaction Matrix normalized by U=",e_t%U_input(ictypat) , &
639           "[Ha] for atom type", e_t%index_typat(ictypat)
640           do iflavor1 = 1, ndim
641             write(e_t%ofile,'(14(a21,2x))',advance="no") (/ ( "...", iflavor2=1,iflavor1 ) /)
642             do iflavor2 = iflavor1+1, ndim
643               icouple = icouple + 1
644               write(e_t%ofile,'(14(es21.14,2x))',advance="no") e_t%uij(icouple,ictypat)
645             end do
646             write(e_t%ofile,*)
647           end do
648         end do
649       end if
650     end do
651     close(e_t%ofile)
652   end subroutine entropyDMFT_dump

ABINIT/m_entropyDMFT/entropyDMFT_init [ Functions ]

[ Top ] [ Functions ]

NAME

  entropyDMFT_init

FUNCTION

  FIXME: add description.

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (J. Bieder)
  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

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

      m_scfcv

CHILDREN

SOURCE

148 subroutine entropyDMFT_init(e_t,dt,pawtab,spacecomm,ifilename,ofilename)
149 
150 !Arguments ------------------------------------
151 
152 !This section has been created automatically by the script Abilint (TD).
153 !Do not modify the following lines by hand.
154 #undef ABI_FUNC
155 #define ABI_FUNC 'entropyDMFT_init'
156 !End of the abilint section
157 
158     type(entropyDMFT_t) , intent(inout) :: e_t
159     type(dataset_type)  , intent(in   ) :: dt
160     type(pawtab_type)   , intent(in   ) :: pawtab(:)
161     integer             , intent(in   ) :: spacecomm
162     character(len=fnlen), intent(in   ) :: ifilename
163     character(len=fnlen), intent(in   ) :: ofilename
164 !Local variables ------------------------------
165     logical :: doRestart
166     integer :: natom, ncatom
167     integer :: iatom, icatom
168     integer :: nctypat,itypat, ictypat
169     integer :: ilambda
170     integer, allocatable :: maptypat(:)
171     character(len=500) :: message
172 
173     e_t%action = dt%dmft_entropy
174     e_t%mylambda = e_t%action-1   ! should be usefull to start form a given value of lambda
175                                   ! -1 added to start from 1 at the first call
176                                   ! of nextLambda
177 
178     if ( e_t%action == AC_NOTHING ) then
179       e_t%isset = .TRUE.
180       return
181     endif
182 
183     e_t%action = AC_ETOT
184 
185     natom = dt%natom
186     ncatom = dt%natpawu
187 
188     icatom = 0
189     do iatom=1,natom
190       if ( pawtab(dt%typat(iatom))%lpawu /= -1 ) icatom = icatom + 1
191     end do
192 
193     if ( icatom /= ncatom ) MSG_ERROR("Inconsistent number of correlated atoms")
194 
195     nctypat = 0
196     do itypat=1,dt%ntypat
197       if ( pawtab(itypat)%lpawu /= -1 ) nctypat = nctypat + 1
198     end do
199 
200     e_t%natom  = natom
201     e_t%ncatom = ncatom
202     e_t%ntypat = dt%ntypat
203     e_t%nctypat = nctypat
204 
205     if ( dt%dmft_nlambda < 3 ) then
206       write(message,'(2a,i4,2a)') "DMFT must have dmft_nlamda >= 3 to compute entropy", &
207         "whereas its value is dmft_nlambda = ",dt%dmft_nlambda,ch10,&
208         "Action : check you input variable dmft_nlambda"
209       MSG_ERROR(message)
210     end if
211 
212     e_t%nlambda = dt%dmft_nlambda
213 
214     doRestart = .false.
215     if ( e_t%nlambda < (e_t%mylambda+1) ) then
216       MSG_ERROR("Restart calculation of DMFT entropy with a value of dmft_entropy greater than dmft_nlambda")
217     else if ( e_t%mylambda > 0 ) then
218       doRestart = .true.
219     end if
220 
221     call entropyDMFT_allocateAll(e_t)
222 
223     e_t%entropy0      = zero
224     e_t%energies(:,:) = zero
225 
226     e_t%lambda(:)     = (/ (DBLE(ilambda-1)/DBLE(e_t%nlambda-1), ilambda=1,e_t%nlambda) /)
227 
228     e_t%temp          = dt%tsmear
229 
230     ! Save each value of U and J for each correlated type
231     ictypat = 1
232     ABI_ALLOCATE(maptypat,(1:dt%ntypat))
233     do itypat=1,dt%ntypat
234       if ( pawtab(itypat)%lpawu /= -1 ) then
235         e_t%lpawu(ictypat) = pawtab(itypat)%lpawu
236         if(dt%dmft_t2g==1.and.e_t%lpawu(ictypat)==2) then
237           e_t%lpawu(ictypat)=1
238         end if
239         e_t%U_input(ictypat) = pawtab(itypat)%upawu
240         e_t%J_input(ictypat) = pawtab(itypat)%jpawu
241         e_t%index_typat(ictypat) = itypat
242         maptypat(itypat) = ictypat
243         ictypat = ictypat + 1
244       end if
245     end do
246 
247     ! Save type local and global of each correlated atom
248     ! Get the correct value of lpawu for each type
249     icatom = 1
250     do iatom=1,e_t%ncatom
251       if ( pawtab(dt%typat(iatom))%lpawu /= -1 ) then
252         e_t%typat(icatom) = maptypat(dt%typat(iatom))
253         e_t%index_atom(icatom) = iatom
254         icatom = icatom + 1
255       end if
256     end do
257     ABI_DEALLOCATE(maptypat)
258 
259     write(message,'(a,1x,78a)') ch10,"+",(/ ("-",ilambda=1,76) /), "+"
260     call wrtout(std_out,message,"COLL")
261     call wrtout(ab_out,message,"COLL")
262     write(message,'(1x,a)') "|             Calculation of entropy within  the DMFT Framework              |"
263     call wrtout(std_out,message,"COLL")
264     call wrtout(ab_out,message,"COLL")
265     write(message,'(1x,40a)') "+",(/ ("- ",ilambda=1,38) /), "+"
266     call wrtout(std_out,message,"COLL")
267     do ilambda = 1, e_t%nlambda
268       write(message,'(1x,a,i4,a11,f6.4,55x,a)') &
269         "|", ilambda, ") lambda = ", e_t%lambda(ilambda), "|"
270       call wrtout(std_out,message,"COLL")
271       do ictypat=1,e_t%nctypat
272         write(message,'(1x,a,6x,a12,i4,4x,a3,5x,2(3x,a4,f6.4,1x,a2),10x,a)') "|", &
273           "- Atom type ", e_t%index_typat(ictypat), "->", &
274           "U = ",e_t%U_input(ictypat)*e_t%lambda(ilambda), "Ha", &
275           "J = ",e_t%J_input(ictypat)*e_t%lambda(ilambda), "Ha","|"
276         call wrtout(std_out,message,"COLL")
277       end do
278     end do
279     write(message,'(1x,78a)') "+",(/ ("-",ilambda=1,76) /), "+"
280     call wrtout(std_out,message,"COLL")
281     call wrtout(ab_out,message,"COLL")
282 
283     ! Set up MPI
284     e_t%spacecomm = spacecomm
285     e_t%rank      = xmpi_comm_rank(spacecomm)
286     e_t%comm_size = xmpi_comm_size(spacecomm)
287 
288     e_t%ofile = get_unit()
289     e_t%ofilename = ofilename
290     e_t%ifilename = ifilename
291     e_t%filename = TRIM(e_t%ofilename)//"_EntropyDMFT"
292     e_t%filedata = TRIM(e_t%ofilename)//"_data4EntropyDMFT"
293 
294     if ( doRestart .eqv. .true. ) then
295       ! If restart fails, then nothing changes and the full calculation is
296       ! perform. Otherwise, we complete as much a possible the structure.
297       call entropyDMFT_restart(e_t)
298     end if
299 
300     ! Rewrite the files with the previous data
301     call entropyDMFT_dump(e_t)
302 
303     e_t%isset = .TRUE.
304 
305 
306   end subroutine entropyDMFT_init

ABINIT/m_entropyDMFT/entropyDMFT_integrate [ Functions ]

[ Top ] [ Functions ]

NAME

  entropyDMFT_integrate

FUNCTION

  FIXME: add description.

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (J. Bieder)
  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

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

      m_entropyDMFT

CHILDREN

SOURCE

1110   subroutine entropyDMFT_integrate(e_t,integrand,integral)
1111 
1112 !Arguments ------------------------------------
1113 
1114 !This section has been created automatically by the script Abilint (TD).
1115 !Do not modify the following lines by hand.
1116 #undef ABI_FUNC
1117 #define ABI_FUNC 'entropyDMFT_integrate'
1118 !End of the abilint section
1119 
1120     type(entropyDMFT_t), intent(inout) :: e_t
1121     real(dp)        , intent(in   ) :: integrand(:,:)
1122     real(dp)        , intent(  out) :: integral
1123 !Local variables ------------------------------
1124     real(dp), allocatable :: ypp(:)
1125     real(dp), allocatable :: fitx(:)
1126     real(dp), allocatable :: fity(:)
1127     real(dp) :: integral1
1128     !real(dp) :: integral2
1129     !real(dp) :: dx
1130     integer :: unit
1131     integer :: Nfit
1132     integer :: icatom
1133     integer :: i
1134 
1135     Nfit = 1000
1136 
1137     ABI_ALLOCATE(ypp,(1:e_t%nlambda))
1138     ABI_ALLOCATE(fitx,(1:Nfit))
1139     ABI_ALLOCATE(fity,(1:Nfit))
1140 
1141     integral = zero
1142 
1143     unit = get_unit()
1144     if ( e_t%rank == 0 ) open(unit=unit,file="fit.log")
1145     do icatom = 1, e_t%ncatom
1146       integral1 = zero
1147       !integral2 = zero
1148       !dx = e_t%U_input(icatom)/dble(e_t%nlambda-1)
1149 
1150       call spline(e_t%lambda,integrand(:,icatom),e_t%nlambda,&
1151         (integrand(2,icatom)-integrand(1,icatom))/e_t%lambda(2),&   ! first derivative left
1152         (integrand(e_t%nlambda,icatom)-integrand(e_t%nlambda-1,icatom))/e_t%lambda(2),& ! first derivative right
1153         ypp)
1154       fitx(1:Nfit)=(/ (dble(i-1)/dble(Nfit-1),i=1,Nfit) /)
1155       call splint(e_t%nlambda,e_t%lambda,integrand(:,icatom),ypp,Nfit,fitx,fity)
1156       integral1 = (sum(fity(1:Nfit))-(fity(1)+fity(Nfit))*half)*e_t%U_input(e_t%typat(icatom))/dble(Nfit-1)
1157       if ( e_t%rank == 0 ) then
1158         do i=1,Nfit
1159           write(unit,'(2ES21.14)') fitx(i),fity(i)
1160         end do
1161       end if
1162       !Unfortunately I don't trust this function and I don't have time to understand it.
1163       !call spline_integrate(integral2,Nfit,dx,integrand(:,icatom))
1164       !write(*,*) integral1, integral2
1165 
1166       !if ( abs(integral2-integral1) >= tol6 ) then
1167       !  write(msg,'(1x,a,1x,f8.6,1x,a,1x,f8.6,1x,a,i4)') "Difference between two different ways of integration is", abs(integral2-integral1), &
1168       !    "which is greater than", tol6, "for correlated atom",e_t%index_atom(icatom)
1169       !  MSG_WARNING(msg)
1170       !end if
1171 
1172       integral = integral+integral1
1173     end do
1174     if ( e_t%rank == 0 ) close(unit)
1175 
1176     ABI_DEALLOCATE(ypp)
1177     ABI_DEALLOCATE(fitx)
1178     ABI_DEALLOCATE(fity)
1179 
1180   end subroutine entropyDMFT_integrate

ABINIT/m_entropyDMFT/entropyDMFT_nextLambda [ Functions ]

[ Top ] [ Functions ]

NAME

  entropyDMFT_nextLambda

FUNCTION

  FIXME: add description.

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (J. Bieder)
  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

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

686   function entropyDMFT_nextLambda(e_t,dt,pawtab,pawang,pawrad) result(nextstep)
687 
688 !Arguments ------------------------------------
689 
690 !This section has been created automatically by the script Abilint (TD).
691 !Do not modify the following lines by hand.
692 #undef ABI_FUNC
693 #define ABI_FUNC 'entropyDMFT_nextLambda'
694 !End of the abilint section
695 
696     type(entropyDMFT_t) , intent(inout) :: e_t
697     type(dataset_type)  , intent(in) :: dt
698     type(pawtab_type)   , intent(inout) :: pawtab(:)
699     type(pawang_type)   , intent(in   ) :: pawang
700     type(pawrad_type)   , intent(inout) :: pawrad(:)
701 !Local variables ------------------------------
702     logical :: nextstep
703     integer :: itypat
704     integer :: mylambda
705     real(dp),allocatable :: upawu(:),jpawu(:)
706     character(len=100) :: message
707 
708     if ( e_t%isset .eqv. .FALSE. ) &
709       MSG_ERROR("entropyDMFT is not initialized")
710 
711     ! go to next lambda
712     mylambda = e_t%mylambda + 1
713     !if ( present(ilambda) ) mylambda = ilambda
714     e_t%mylambda = mylambda
715 
716     if ( e_t%action == AC_NOTHING .and. mylambda >= 1 ) then
717       nextstep = .FALSE.
718       ! We do nothing and return false, one scfvc has already performed
719     else if ( e_t%action == AC_NOTHING .and. mylambda < 1 ) then
720       nextstep = .TRUE.
721       ! We do nothing but return true to perform scfcv a least once
722     else if ( e_t%action == AC_ETOT .and. mylambda <= e_t%nlambda ) then ! we iterate over lambda
723       nextstep = .TRUE.
724       ABI_ALLOCATE(upawu,(dt%ntypat))
725       ABI_ALLOCATE(jpawu,(dt%ntypat))
726       do itypat = 1, e_t%nctypat
727         upawu(e_t%index_typat(itypat)) = e_t%lambda(mylambda) * e_t%U_input(itypat)
728         jpawu(e_t%index_typat(itypat)) = e_t%lambda(mylambda) * e_t%J_input(itypat)
729       end do
730     else ! we did all lambda values
731       nextstep = .FALSE.
732     endif
733 
734     if ( e_t%action == AC_ETOT .and. nextstep .eqv. .true. ) then
735       write(message,'(a,1x,78a)') ch10,"+",(/ ("-",itypat=1,76) /), "+"
736       call wrtout(std_out,message,"COLL")
737       call wrtout(ab_out,message,"COLL")
738       write(message,'(1x,a,i4,a11,f6.4,55x,a)') &
739         "|", mylambda, ") lambda = ", e_t%lambda(mylambda), "|"
740       call wrtout(std_out,message,"COLL")
741       call wrtout(ab_out,message,"COLL")
742       do itypat=1,e_t%nctypat
743         write(message,'(1x,a,6x,a12,i4,4x,a3,5x,2(3x,a4,f6.4,1x,a2),10x,a)') "|",&
744           "- Atom type ", e_t%index_typat(itypat), "->", &
745           "U = ",e_t%U_input(itypat)*e_t%lambda(mylambda), "Ha", &
746           "J = ",e_t%J_input(itypat)*e_t%lambda(mylambda), "Ha","|"
747           call wrtout(std_out,message,"COLL")
748           call wrtout(ab_out,message,"COLL")
749       end do
750       write(message,'(1x,78a)') "+",(/ ("-",itypat=1,76) /), "+"
751       call wrtout(std_out,message,"COLL")
752       call wrtout(ab_out,message,"COLL")
753       call pawpuxinit(dt%dmatpuopt,dt%exchmix,dt%f4of2_sla,dt%f6of2_sla,&
754 &        jpawu,dt%lexexch,dt%lpawu,dt%ntypat,pawang,dt%pawprtvol,&
755 &        pawrad,pawtab,upawu,dt%usedmft,dt%useexexch,dt%usepawu)
756       ABI_DEALLOCATE(upawu)
757       ABI_DEALLOCATE(jpawu)
758     end if
759 
760   end function entropyDMFT_nextLambda

ABINIT/m_entropyDMFT/entropyDMFT_restart [ Functions ]

[ Top ] [ Functions ]

NAME

  entropyDMFT_restart

FUNCTION

  FIXME: add description.

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (J. Bieder)
  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

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

      m_entropyDMFT

CHILDREN

SOURCE

404   subroutine entropyDMFT_restart(e_t)
405 
406 !Arguments ------------------------------------
407 
408 !This section has been created automatically by the script Abilint (TD).
409 !Do not modify the following lines by hand.
410 #undef ABI_FUNC
411 #define ABI_FUNC 'entropyDMFT_restart'
412 !End of the abilint section
413 
414     type(entropyDMFT_t), intent(inout) :: e_t
415 !Local variables-------------------------------
416     logical            :: doBcast
417     character(len=200) :: msg
418     character(len=21 ) :: hdr
419     integer            :: iostate
420     integer            :: ilambda
421     integer            :: tlambda
422     integer            :: ictypat
423     integer            :: ndim
424     integer            :: icouple
425     integer            :: iflavor1
426     integer            :: iflavor2
427     integer            :: icatom
428     integer            :: iatom
429     real(dp)           :: lambda
430 
431     write(msg,'(a,1x,2a,1x,a,f6.4,a)') ch10,"EtotDMFT will try to restart calculation from a previous run", ch10, &
432     "and this calculation should start at lambda = ",e_t%lambda(e_t%mylambda+1), ch10
433     call wrtout(std_out,msg,"COLL")
434     call wrtout(ab_out,msg,"COLL")
435 
436     doBcast = .true.
437 
438     if ( e_t%rank == 0 ) then
439       tlambda = 1
440       inquire(file=e_t%filedata,iostat=iostate)
441 
442       if ( iostate /= 0 ) then
443         write(msg,'(5a)') "File ", trim(e_t%filedata), " does not exist or is not accessible", ch10, &
444         "-> No restart performed but full calculation."
445         MSG_WARNING(msg)
446         e_t%mylambda = 0
447       end if
448 
449       open(unit=e_t%ofile,file=e_t%filedata,action="read",form="unformatted")
450       read(e_t%ofile,end=42) hdr, iostate
451 
452       if ( hdr /= HDR_NAME .or. iostate /= 1 ) then
453         write(msg,'(5a)') "File ", trim(e_t%filedata), " does not contain the proper header", ch10, &
454         "-> No restart performed but full calculation."
455         MSG_WARNING(msg)
456         e_t%mylambda = 0
457       end if
458 
459       do ilambda = 1, e_t%mylambda
460         tlambda = ilambda
461         read(e_t%ofile,end=42) lambda
462         if ( ABS(lambda - e_t%lambda(ilambda)) >= tol9 ) then
463           write(msg,'(5a,f6.4,a,f6.4)') "File ", trim(e_t%filedata), " is wrong:", ch10, &
464           "Lambda values are differente: in file ", lambda, " instead of ", e_t%lambda(ilambda)
465           MSG_WARNING(msg)
466           goto 42
467         end if
468 
469         if ( ilambda == 1 ) then
470           read(e_t%ofile,end=42) lambda
471           if ( lambda /= e_t%temp ) then
472             write(msg,'(7a,f6.4)') "File ", trim(e_t%filedata), " is wrong:", ch10, &
473             "Temperature is different than the value of tsmear", ch10, &
474             "-> No restart performed but full calculation."
475             MSG_WARNING(msg)
476             goto 42
477           end if
478           read(e_t%ofile,end=42) e_t%entropy0
479           read(e_t%ofile,end=42) e_t%energies(E_DC,E_U0)
480         else if ( ilambda == e_t%nlambda ) then ! should never happend ?!
481           read(e_t%ofile,end=42) e_t%energies(E_DC,E_UU)
482           do ictypat = 1, e_t%nctypat
483             ndim = 2*(2*e_t%lpawu(ictypat)+1)
484             icouple = 0
485             do iflavor1 = 1, ndim
486               do iflavor2 = iflavor1+1, ndim
487                 icouple = icouple + 1
488                 read(e_t%ofile,end=42) e_t%uij(icouple,ictypat)
489               end do
490             end do
491           end do
492         end if
493 
494         do icatom = 1, e_t%ncatom
495           iatom = e_t%index_atom(icatom)
496           ndim = 2*(2*e_t%lpawu(e_t%typat(icatom))+1)
497           icouple = 0
498           read(e_t%ofile,end=42) e_t%e_dc(icatom,ilambda)
499           do iflavor1 = 1, ndim
500             do iflavor2 = iflavor1+1, ndim
501               icouple = icouple + 1
502               read(e_t%ofile,end=42) e_t%docc(icouple,icatom,ilambda)
503             end do
504           end do
505         end do
506       end do
507       close(e_t%ofile)
508       goto 43
509 42    write(msg,'(5a,f6.4)') "File ", trim(e_t%filedata), " is wrong or incomplete", ch10, &
510           "-> Restart calculation will restart at lambda = ",e_t%lambda(tlambda)
511       MSG_WARNING(msg)
512       close(e_t%ofile)
513       e_t%mylambda = tlambda-1 ! -1 to go to previous lambda
514     end if
515     ! MPI BDCAST
516 43  call xmpi_bcast(e_t%mylambda,0, e_t%spacecomm, ictypat)
517     call xmpi_bcast(e_t%entropy0,0, e_t%spacecomm, ictypat)
518     call xmpi_bcast(e_t%energies,0, e_t%spacecomm, ictypat)
519     !call xmpi_bcast(e_t%uij,0, e_t%spacecomm, ictypat) ! No need since it
520     ! restart always perform the last lambda and this is calculation at the end
521     ! of last lambda
522     call xmpi_bcast(e_t%e_dc,0, e_t%spacecomm, ictypat)
523     call xmpi_bcast(e_t%docc,0, e_t%spacecomm, ictypat)
524 
525   end subroutine entropyDMFT_restart

m_entropyDMFT/entropyDMFT [ Types ]

[ Top ] [ m_entropyDMFT ] [ Types ]

NAME

  entropyDMFT

FUNCTION

  This structured datatype contains the necessary data

COPYRIGHT

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

SOURCE

 81   type, public :: entropyDMFT_t
 82     logical               :: isset=.FALSE.  ! flag to be sure we are initialize
 83     integer               :: spacecomm      ! MPI comm
 84     integer               :: rank           ! rank in the comm
 85     integer               :: comm_size      ! Number of cpus in the comm
 86     integer               :: action         ! what to do in gstate
 87     integer               :: mylambda       ! Current lambda
 88     integer               :: natom          ! number of atoms
 89     integer               :: ncatom         ! number of correlated atoms
 90     integer               :: ntypat         ! number of type of atoms
 91     integer               :: nctypat        ! number of type of correlated atoms
 92     integer               :: nlambda        ! number of integration points
 93     integer               :: ofile          ! unit for file output data
 94     character(len=fnlen)  :: filename       ! name fo the file user readable
 95     character(len=fnlen)  :: filedata       ! name fo the file for restart purposes
 96     character(len=fnlen)  :: ofilename      ! ofstream prefix
 97     character(len=fnlen)  :: ifilename      ! ifstream prefix
 98     real(dp)              :: temp           ! temperature
 99     real(dp)              :: entropy0       ! entropy for lambda=0
100     real(dp)              :: energies(2,2)  ! internal energy for lambda=0 and 1
101     integer , allocatable :: index_atom(:)  ! index for correlated atoms
102     integer , allocatable :: index_typat(:) ! index for correlated types
103     integer , allocatable :: lpawu(:)       ! orbital moment to treat (ntypat,2)
104     integer , allocatable :: typat(:)       ! type of each correlated atom
105     real(dp), allocatable :: U_input(:)     ! U from input file (ntypat)
106     real(dp), allocatable :: J_input(:)     ! J from input file (ntypat)
107     real(dp), allocatable :: lambda(:)      ! ilamda
108     real(dp), allocatable :: docc(:,:,:)    ! n_In_j, natom, ilamda
109     real(dp), allocatable :: e_dc(:,:)      ! natom, ilamda
110     real(dp), allocatable :: uij(:,:)       ! uij to compute <uij n_i n_j>, ntypat
111   end type entropyDMFT_t