TABLE OF CONTENTS


ABINIT/m_sigma [ Modules ]

[ Top ] [ Modules ]

NAME

  m_sigma

FUNCTION

  This module provides the definition of the sigma_t data type
  used to store results of the GW calculation as well as as
  methods bound to the object.

COPYRIGHT

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

PARENTS

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 MODULE m_sigma
28 
29  use defs_basis
30  use defs_datatypes
31  use defs_abitypes
32  use m_xmpi
33  use m_abicore
34  use m_errors
35  use iso_c_binding
36  use m_nctk
37 #ifdef HAVE_NETCDF
38  use netcdf
39 #endif
40  use m_wfd
41 
42  use m_numeric_tools,  only : c2r
43  use m_gwdefs,         only : unt_gw, unt_sig, unt_sgr, unt_sgm, unt_gwdiag, sigparams_t, sigma_needs_w
44  use m_crystal,        only : crystal_t
45  use m_bz_mesh,        only : kmesh_t, littlegroup_t, findqg0
46  use m_screening,      only : epsilonm1_results
47 
48  implicit none
49 
50  private

m_sigma/find_wpoles_for_cd [ Functions ]

[ Top ] [ m_sigma ] [ Functions ]

NAME

  find_wpoles_for_cd

FUNCTION

  Find the max frequency needed to account for all the poles
  of GW used in the contour deformation technique.

INPUTS

  Sigp=sigparams_t

OUTPUT

  omega_max

NOTES

PARENTS

CHILDREN

      findqg0,wfd_distribute_bands,wfd_update_bkstab,xmpi_sum

SOURCE

1301 subroutine find_wpoles_for_cd(Sigp,Sr,Kmesh,BSt,omega_max)
1302 
1303 
1304 !This section has been created automatically by the script Abilint (TD).
1305 !Do not modify the following lines by hand.
1306 #undef ABI_FUNC
1307 #define ABI_FUNC 'find_wpoles_for_cd'
1308 !End of the abilint section
1309 
1310  implicit none
1311 
1312 !Arguments ------------------------------------
1313 !scalars
1314  type(sigparams_t),intent(in) :: Sigp
1315  type(sigma_t),intent(in) :: Sr
1316  type(ebands_t),intent(in) :: Bst
1317  type(kmesh_t),intent(in) :: Kmesh
1318  real(dp),intent(out) :: omega_max
1319 
1320 !Local variables-------------------------------
1321 !scalars
1322  integer :: spin,ik_ibz,band_gr,bgw_start,bgw_stop,io,ioe0j
1323  integer :: ikgw,ikgw_ibz,ikgw_bz,band_gw,nomega_tot
1324  real(dp) :: e_green,e_screen,theta_mu_minus_e0i,e_qp
1325  real(dp) :: fact_sp
1326  !character(len=500) :: msg
1327 !arrays
1328  real(dp),allocatable :: omegame0i(:)
1329 
1330 ! *************************************************************************
1331 
1332  omega_max = smallest_real
1333  !
1334  ! === Normalization of theta_mu_minus_e0i ===
1335  ! * If nsppol==2, qp_occ $\in [0,1]$
1336  fact_sp=one
1337  if (Bst%nsppol==1) then
1338    fact_sp=half; if (Bst%nspinor==2) fact_sp=one
1339  end if
1340  !
1341  ! Total number of frequencies for sigma (Spectral function + mesh for the derivative).
1342  nomega_tot=Sr%nomega_r+Sr%nomega4sd
1343  ABI_MALLOC(omegame0i,(nomega_tot))
1344 
1345  ioe0j=Sr%nomega4sd/2+1
1346  !
1347  ! Loop over bands used to construct the Green function.
1348  do spin=1,Bst%nsppol
1349    do ik_ibz=1,Bst%nkpt
1350      do band_gr=1,Bst%nband(ik_ibz+(spin-1)*Bst%nkpt)
1351        e_green           = Bst%eig(band_gr,ik_ibz,spin)
1352        theta_mu_minus_e0i= Bst%occ(band_gr,ik_ibz,spin)*fact_sp
1353        !
1354        ! Loop over GW states.
1355        do ikgw=1,Sigp%nkptgw
1356          bgw_start=Sigp%minbnd(ikgw,spin)
1357          bgw_stop =Sigp%minbnd(ikgw,spin)
1358          ikgw_bz  =Sigp%kptgw2bz(ikgw_bz)
1359          ikgw_ibz =Kmesh%tab(ikgw_bz)
1360 
1361          do band_gw=bgw_start,bgw_stop
1362            e_qp      = Bst%eig(band_gw,ikgw_ibz,spin)
1363            !
1364            ! Get frequencies $\omega$-\epsilon_in$ to evaluate  $d\Sigma/dE$, note the spin
1365            ! subtract e_KS since we have stored e_KS+ Delta \omega in Sr%omega4sd, not required for AC
1366            if (Sr%nomega_r>0) omegame0i(1:Sr%nomega_r)=DBLE(Sigp%omega_r(1:Sr%nomega_r))-e_green
1367            do io=Sr%nomega_r+1,nomega_tot
1368              !omegame0i(io)=DBLE(Sr%omega4sd(band_gw,ikgw_ibz,io-Sr%nomega_r,spin)) - e_green
1369              !Sr%omega4sd(jb,ik_ibz,io,spin)=Sr%egw(jb,ik_ibz,spin)+Sigp%deltae*(io-ioe0j)
1370              omegame0i(io) = e_qp + Sigp%deltae*(io-ioe0j) - e_green
1371            end do
1372 
1373            do io=1,nomega_tot
1374              e_screen =  ABS(omegame0i(io))
1375              if (omegame0i(io)>tol12) then
1376                !ket(spadc+ig,ios)=ket(spadc+ig,ios)+ct*(one-theta_mu_minus_e0i)
1377                if ( (one-theta_mu_minus_e0i) > tol12 ) omega_max = MAX(omega_max, e_screen)
1378              end if
1379              if (omegame0i(io)<-tol12) then
1380                !ket(spadc+ig,ios)=ket(spadc+ig,ios)-ct*theta_mu_minus_e0i
1381                if ( theta_mu_minus_e0i > tol12) omega_max = MAX(omega_max, e_screen)
1382              end if
1383            end do
1384 
1385          end do
1386        end do
1387        !
1388      end do
1389    end do
1390  end do
1391 
1392  ABI_FREE(omegame0i)
1393 
1394 end subroutine find_wpoles_for_cd

m_sigma/gw_spectral_function [ Functions ]

[ Top ] [ m_sigma ] [ Functions ]

NAME

 gw_spectral_function

FUNCTION

  Compute the spectral function in the GW approximation

INPUTS

  io,ib,ikibz,is=Frequency, band, k-point, spin index
  Sr=sigma results datatype

OUTPUT

PARENTS

SOURCE

622 function gw_spectral_function(Sr,io,ib,ikibz,is) result(aw)
623 
624 
625 !This section has been created automatically by the script Abilint (TD).
626 !Do not modify the following lines by hand.
627 #undef ABI_FUNC
628 #define ABI_FUNC 'gw_spectral_function'
629 !End of the abilint section
630 
631  implicit none
632 
633 !Arguments ------------------------------------
634 !scalars
635  integer,intent(in) :: io,ib,ikibz,is
636  real(dp) :: aw
637  type(sigma_t),intent(in) :: Sr
638 
639 ! *********************************************************************
640 
641  aw = one/pi*ABS(AIMAG(Sr%sigcme(ib,ikibz,io,is)))&
642 &  /( (REAL(Sr%omega_r(io)-Sr%hhartree(ib,ib,ikibz,is)-Sr%sigxcme(ib,ikibz,io,is)))**2&
643 &    +(AIMAG(Sr%sigcme(ib,ikibz,io,is)))**2) /Ha_eV
644 
645 end function gw_spectral_function

m_sigma/print_Sigma_perturbative [ Functions ]

[ Top ] [ m_sigma ] [ Functions ]

NAME

 print_Sigma_perturbative

FUNCTION

  Write the results of the GW calculation done with the perturbative approach

INPUTS

OUTPUT

PARENTS

      m_sigma

CHILDREN

      findqg0,wfd_distribute_bands,wfd_update_bkstab,xmpi_sum

SOURCE

669 subroutine print_Sigma_perturbative(Sr,ik_ibz,iband,isp,unit,prtvol,mode_paral,witheader)
670 
671 
672 !This section has been created automatically by the script Abilint (TD).
673 !Do not modify the following lines by hand.
674 #undef ABI_FUNC
675 #define ABI_FUNC 'print_Sigma_perturbative'
676 !End of the abilint section
677 
678  implicit none
679 
680 !Arguments ------------------------------------
681 !scalars
682  integer,intent(in) :: iband,ik_ibz,isp
683  integer,optional,intent(in) :: prtvol,unit
684  character(len=*),optional,intent(in) :: mode_paral
685  logical,optional,intent(in) :: witheader
686  type(sigma_t),intent(in) :: Sr
687 
688 !Local variables-------------------------------
689 !scalars
690  integer :: my_unt,verbose
691  character(len=4) :: my_mode
692  character(len=500) :: msg
693 
694 ! *********************************************************************
695 
696  my_unt =std_out; if (PRESENT(unit      )) my_unt =unit
697  verbose=0      ; if (PRESENT(prtvol    )) verbose=prtvol
698  my_mode='COLL' ; if (PRESENT(mode_paral)) my_mode=mode_paral
699 
700  if (PRESENT(witheader)) then
701    if (witheader) then
702      call wrtout(my_unt,' Band     E0 <VxcLDA>   SigX SigC(E0)      Z dSigC/dE  Sig(E)    E-E0       E ',my_mode)
703    end if
704  end if
705 
706  if (Sr%usepawu==0) then
707 
708    if (Sr%nsig_ab/=1) then
709      write(msg,'(i5,9f8.3)')                        &
710 &           iband,                                  &
711 &           Sr%e0          (iband,ik_ibz,1)*Ha_eV,  &
712 &           SUM(Sr%vxcme   (iband,ik_ibz,:))*Ha_eV, &
713 &           SUM(Sr%sigxme  (iband,ik_ibz,:))*Ha_eV, &
714 &      REAL(SUM(Sr%sigcmee0(iband,ik_ibz,:)))*Ha_eV,&
715 &      REAL(Sr%ze0         (iband,ik_ibz,1)),       &
716 &      REAL(SUM(Sr%dsigmee0(iband,ik_ibz,:))),      &
717 &      REAL(SUM(Sr%sigmee  (iband,ik_ibz,:)))*Ha_eV,&
718 &      REAL(Sr%degw        (iband,ik_ibz,1))*Ha_eV, &
719 &      REAL(Sr%egw         (iband,ik_ibz,1))*Ha_eV
720        call wrtout(my_unt,msg,my_mode)
721      if (verbose/=0) then
722        write(msg,'(i5,9f8.3)')                         &
723 &              iband,                                  &
724 &              zero,                                   &
725 &              zero,                                   &
726 &              zero,                                   &
727 &        AIMAG(SUM(Sr%sigcmee0(iband,ik_ibz,:)))*Ha_eV,&
728 &        AIMAG(Sr%ze0         (iband,ik_ibz,1)),       &
729 &        AIMAG(SUM(Sr%dsigmee0(iband,ik_ibz,:))),      &
730 &        AIMAG(SUM(Sr%sigmee  (iband,ik_ibz,:)))*Ha_eV,&
731 &        AIMAG(Sr%degw        (iband,ik_ibz,1))*Ha_eV, &
732 &        AIMAG(Sr%egw         (iband,ik_ibz,1))*Ha_eV
733        call wrtout(my_unt,msg,my_mode)
734      end if
735   else
736     write(msg,'(i5,9f8.3)')                     &
737 &          iband,                               &
738 &          Sr%e0      (iband,ik_ibz,isp)*Ha_eV, &
739 &          Sr%vxcme   (iband,ik_ibz,isp)*Ha_eV, &
740 &          Sr%sigxme  (iband,ik_ibz,isp)*Ha_eV, &
741 &     REAL(Sr%sigcmee0(iband,ik_ibz,isp))*Ha_eV,&
742 &     REAL(Sr%ze0     (iband,ik_ibz,isp)),      &
743 &     REAL(Sr%dsigmee0(iband,ik_ibz,isp)),      &
744 &     REAL(Sr%sigmee  (iband,ik_ibz,isp))*Ha_eV,&
745 &     REAL(Sr%degw    (iband,ik_ibz,isp))*Ha_eV,&
746 &     REAL(Sr%egw     (iband,ik_ibz,isp))*Ha_eV
747     call wrtout(my_unt,msg,my_mode)
748 
749     if (verbose/=0) then
750       write(msg,'(i5,9f8.3)')                       &
751 &              iband,                               &
752 &              zero,                                &
753 &              zero,                                &
754 &              zero,                                &
755 &        AIMAG(Sr%sigcmee0(iband,ik_ibz,isp))*Ha_eV,&
756 &        AIMAG(Sr%ze0     (iband,ik_ibz,isp)),      &
757 &        AIMAG(Sr%dsigmee0(iband,ik_ibz,isp)),      &
758 &        AIMAG(Sr%sigmee  (iband,ik_ibz,isp))*Ha_eV,&
759 &        AIMAG(Sr%degw    (iband,ik_ibz,isp))*Ha_eV,&
760 &        AIMAG(Sr%egw     (iband,ik_ibz,isp))*Ha_eV
761        call wrtout(my_unt,msg,my_mode)
762     end if
763   end if
764 
765  else  ! PAW+U+GW calculation.
766    ABI_CHECK(Sr%nsig_ab==1,'LDA+U with spinor not implemented')
767    write(msg,'(i5,10f8.3)')                    &
768 &         iband,                               &
769 &         Sr%e0      (iband,ik_ibz,isp)*Ha_eV, &
770 &         Sr%vxcme   (iband,ik_ibz,isp)*Ha_eV, &
771 &         Sr%vUme    (iband,ik_ibz,isp)*Ha_eV, &
772 &         Sr%sigxme  (iband,ik_ibz,isp)*Ha_eV, &
773 &    REAL(Sr%sigcmee0(iband,ik_ibz,isp))*Ha_eV,&
774 &    REAL(Sr%ze0     (iband,ik_ibz,isp)),      &
775 &    REAL(Sr%dsigmee0(iband,ik_ibz,isp)),      &
776 &    REAL(Sr%sigmee  (iband,ik_ibz,isp))*Ha_eV,&
777 &    REAL(Sr%degw    (iband,ik_ibz,isp))*Ha_eV,&
778 &    REAL(Sr%egw     (iband,ik_ibz,isp))*Ha_eV
779    call wrtout(my_unt,msg,my_mode)
780 
781    if (verbose/=0) then
782      write(msg,'(i5,10f8.3)')                    &
783 &           iband,                               &
784 &           zero,                                &
785 &           zero,                                &
786 &           zero,                                &
787 &           zero,                                &
788 &     AIMAG(Sr%sigcmee0(iband,ik_ibz,isp))*Ha_eV,&
789 &     AIMAG(Sr%ze0     (iband,ik_ibz,isp)),      &
790 &     AIMAG(Sr%dsigmee0(iband,ik_ibz,isp)),      &
791 &     AIMAG(Sr%sigmee  (iband,ik_ibz,isp))*Ha_eV,&
792 &     AIMAG(Sr%degw    (iband,ik_ibz,isp))*Ha_eV,&
793 &     AIMAG(Sr%egw     (iband,ik_ibz,isp))*Ha_eV
794       call wrtout(my_unt,msg,my_mode)
795    end if
796  end if
797 
798 end subroutine print_Sigma_perturbative

m_sigma/print_Sigma_QPSC [ Functions ]

[ Top ] [ m_sigma ] [ Functions ]

NAME

  print_Sigma_QPSC

FUNCTION

  Write the results of the GW calculation in case of self-consistency

INPUTS

OUTPUT

PARENTS

      m_sigma

CHILDREN

      findqg0,wfd_distribute_bands,wfd_update_bkstab,xmpi_sum

SOURCE

822 subroutine print_Sigma_QPSC(Sr,ik_ibz,iband,isp,KS_BSt,unit,prtvol,mode_paral)
823 
824 
825 !This section has been created automatically by the script Abilint (TD).
826 !Do not modify the following lines by hand.
827 #undef ABI_FUNC
828 #define ABI_FUNC 'print_Sigma_QPSC'
829 !End of the abilint section
830 
831  implicit none
832 
833 !Arguments ------------------------------------
834 !scalars
835  integer,intent(in) :: iband,ik_ibz,isp
836  integer,intent(in),optional :: prtvol,unit
837  character(len=*),intent(in),optional :: mode_paral
838  type(sigma_t),intent(in) :: Sr
839  type(ebands_t),intent(in) :: KS_BSt
840 
841 !Local variables-------------------------------
842 !scalars
843  integer :: my_unt,verbose
844  character(len=4) :: my_mode
845  character(len=500) :: msg
846 
847 ! *********************************************************************
848 
849  my_unt =std_out; if (PRESENT(unit      )) my_unt =unit
850  verbose=0      ; if (PRESENT(prtvol    )) verbose=prtvol
851  my_mode='COLL' ; if (PRESENT(mode_paral)) my_mode=mode_paral
852 
853 ! write(msg,'(a)')&
854 !&   ' Band     E_lda   <Vxclda>   E(N-1)  <Hhartree>   SigX  SigC[E(N-1)]',&
855 !&   '    Z     dSigC/dE  Sig[E(N)]  DeltaE  E(N)_pert E(N)_diago'
856 
857  if (Sr%usepawu==0 .or. .TRUE.) then
858    if (Sr%nsig_ab/=1) then
859      write(msg,'(i5,12(2x,f8.3))')                        &
860 &           iband,                                        &
861 &           KS_BSt%eig     (iband,ik_ibz,1)*Ha_eV,        &
862 &           SUM(Sr%vxcme   (iband,ik_ibz,:))*Ha_eV,       &
863 &           Sr%e0          (iband,ik_ibz,1)*Ha_eV,        &
864 &      REAL(SUM(Sr%hhartree(iband,iband,ik_ibz,:)))*Ha_eV,&
865 &           SUM(Sr%sigxme  (iband,ik_ibz,:))*Ha_eV,       &
866 &      REAL(SUM(Sr%sigcmee0(iband,ik_ibz,:)))*Ha_eV,      &
867 &      REAL(Sr%ze0         (iband,ik_ibz,1)),             &
868 &      REAL(SUM(Sr%dsigmee0(iband,ik_ibz,:))),            &
869 &      REAL(SUM(Sr%sigmee  (iband,ik_ibz,:)))*Ha_eV,      &
870 &      REAL(Sr%degw        (iband,ik_ibz,1))*Ha_eV,       &
871 &      REAL(Sr%egw         (iband,ik_ibz,1))*Ha_eV,       &
872 &           Sr%en_qp_diago (iband,ik_ibz,1)*Ha_eV
873      call wrtout(my_unt,msg,my_mode)
874 
875      write(msg,'(i5,12(2x,f8.3))')                         &
876 &            iband,                                        &
877 &            zero,                                         &
878 &            zero,                                         &
879 &            zero,                                         &
880 &      AIMAG(SUM(Sr%hhartree(iband,iband,ik_ibz,:)))*Ha_eV,&
881 &            zero,                                         &
882 &      AIMAG(SUM(Sr%sigcmee0(iband,ik_ibz,:)))*Ha_eV,      &
883 &      AIMAG(Sr%ze0         (iband,ik_ibz,1)),             &
884 &      AIMAG(SUM(Sr%dsigmee0(iband,ik_ibz,:))),            &
885 &      AIMAG(SUM(Sr%sigmee  (iband,ik_ibz,:)))*Ha_eV,      &
886 &      AIMAG(Sr%degw        (iband,ik_ibz,1))*Ha_eV,       &
887 &      AIMAG(Sr%egw         (iband,ik_ibz,1))*Ha_eV,       &
888 &            zero
889      if (verbose/=0) then
890        call wrtout(my_unt,msg,my_mode)
891      end if
892    else
893      write(msg,'(i5,12(2x,f8.3))')                        &
894 &           iband,                                        &
895 &           KS_BSt%eig    (iband,ik_ibz,isp)*Ha_eV,       &
896 &           Sr%vxcme      (iband,ik_ibz,isp)*Ha_eV,       &
897 &           Sr%e0         (iband,ik_ibz,isp)*Ha_eV,       &
898 &      REAL(Sr%hhartree   (iband,iband,ik_ibz,isp))*Ha_eV,&
899 &           Sr%sigxme     (iband,ik_ibz,isp)*Ha_eV,       &
900 &      REAL(Sr%sigcmee0   (iband,ik_ibz,isp))*Ha_eV,      &
901 &      REAL(Sr%ze0        (iband,ik_ibz,isp)),            &
902 &      REAL(Sr%dsigmee0   (iband,ik_ibz,isp)),            &
903 &      REAL(Sr%sigmee     (iband,ik_ibz,isp))*Ha_eV,      &
904 &      REAL(Sr%degw       (iband,ik_ibz,isp))*Ha_eV,      &
905 &      REAL(Sr%egw        (iband,ik_ibz,isp))*Ha_eV,      &
906 &           Sr%en_qp_diago(iband,ik_ibz,isp)*Ha_eV
907      call wrtout(my_unt,msg,my_mode)
908 
909      write(msg,'(i5,12(2x,f8.3))')                        &
910 &            iband,                                       &
911 &            zero,                                        &
912 &            zero,                                        &
913 &            zero,                                        &
914 &      AIMAG(Sr%hhartree  (iband,iband,ik_ibz,isp))*Ha_eV,&
915 &            zero,                                        &
916 &      AIMAG(Sr%sigcmee0   (iband,ik_ibz,isp))*Ha_eV,     &
917 &      AIMAG(Sr%ze0        (iband,ik_ibz,isp)),           &
918 &      AIMAG(Sr%dsigmee0   (iband,ik_ibz,isp)),           &
919 &      AIMAG(Sr%sigmee     (iband,ik_ibz,isp))*Ha_eV,     &
920 &      AIMAG(Sr%degw       (iband,ik_ibz,isp))*Ha_eV,     &
921 &      AIMAG(Sr%egw        (iband,ik_ibz,isp))*Ha_eV,     &
922 &            zero
923      if (verbose/=0) then
924        call wrtout(my_unt,msg,my_mode)
925      end if
926    end if
927 
928  else ! PAW+U+GW calculation.
929    MSG_ERROR("PAW+U+GW not yet implemented")
930  end if
931 
932 end subroutine print_Sigma_QPSC

m_sigma/sigma_distribute_bks [ Functions ]

[ Top ] [ m_sigma ] [ Functions ]

NAME

  sigma_distribute_bks

FUNCTION

  Distribute the loop over (b,k,s) used to calculate the self-energy matrix elements
  taking into account the MPI distribution of the wavefunctions and the use of
  symmetries to reduce the BZ sum to an appropriate irreducible wedge.

INPUTS

 nsppol
 can_symmetrize(nsppol)=.TRUE if symmetries can be used to reduce the number of k-points to be summed.
 Kmesh<kmesh_t>
 Qmesh<kmesh_t>
 Ltg_kgw<littlegroup_t>
 Wfd(wfd_t),intent(inout) ::
 mg0(3)
 kptgw(3)
 [bks_mask(Wfd%mband,Kmesh%nbz,nsppol)]
 [got(Wfd%nproc)]=The number of tasks already assigned to the nodes.
 [global]=If true, an MPI global communication is performed such that each node will have the same table. Useful
   if for implementing algorithms in which each node needs to know the global distribution of the tasks, not only
   the task it has to complete. Defaults to .FALSE.

OUTPUT

  my_nbks
  proc_distrb(Wfd%mband,Kmesh%nbz,nsppol)

SIDE EFFECTS

  Wfd%bks_tab

PARENTS

      calc_sigc_me,calc_sigx_me,cohsex_me

CHILDREN

      findqg0,wfd_distribute_bands,wfd_update_bkstab,xmpi_sum

SOURCE

1724 subroutine sigma_distribute_bks(Wfd,Kmesh,Ltg_kgw,Qmesh,nsppol,can_symmetrize,kptgw,mg0,my_nbks,proc_distrb,got,bks_mask,global)
1725 
1726 
1727 !This section has been created automatically by the script Abilint (TD).
1728 !Do not modify the following lines by hand.
1729 #undef ABI_FUNC
1730 #define ABI_FUNC 'sigma_distribute_bks'
1731 !End of the abilint section
1732 
1733  implicit none
1734 
1735 !Arguments ------------------------------------
1736 !scalars
1737  integer,intent(in) :: nsppol
1738  integer,intent(out) :: my_nbks
1739  logical,optional,intent(in) :: global
1740  type(kmesh_t),intent(in) :: Kmesh,Qmesh
1741  type(littlegroup_t),intent(in) :: Ltg_kgw
1742  type(wfd_t),intent(inout) :: Wfd
1743 !arrays
1744  integer,intent(in) :: mg0(3)
1745  integer,optional,intent(inout) :: got(Wfd%nproc)
1746  integer,intent(out) :: proc_distrb(Wfd%mband,Kmesh%nbz,nsppol)
1747  real(dp),intent(in) :: kptgw(3)
1748  logical,intent(in) :: can_symmetrize(Wfd%nsppol)
1749  logical,optional,intent(in) :: bks_mask(Wfd%mband,Kmesh%nbz,nsppol)
1750 
1751 !Local variables-------------------------------
1752 !scalars
1753  integer :: ierr,ik_bz,ik_ibz,spin,iq_bz,my_nband
1754  !character(len=500) :: msg
1755 !arrays
1756  integer :: g0(3)
1757  real(dp) :: kgwmk(3)
1758  integer :: get_more(Wfd%nproc),my_band_list(Wfd%mband)
1759  !integer :: test(Wfd%mband,Kmesh%nbz,nsppol)
1760  logical :: bmask(Wfd%mband)
1761 
1762 !************************************************************************
1763 
1764  call wfd_update_bkstab(Wfd)
1765 
1766  get_more=0; if (PRESENT(got)) get_more=got
1767 
1768  ! Different distribution of tasks depending whether symmetries can be used or not.
1769  proc_distrb= xmpi_undefined_rank
1770 
1771  do spin=1,Wfd%nsppol
1772 
1773    if (can_symmetrize(spin)) then
1774      do ik_bz=1,Kmesh%nbz
1775        ik_ibz = Kmesh%tab(ik_bz)
1776        kgwmk= kptgw-Kmesh%bz(:,ik_bz) ! kptgw must be inside the BZ
1777        call findqg0(iq_bz,g0,kgwmk,Qmesh%nbz,Qmesh%bz,mG0) ! Identify q_bz and G0 where q_bz+G0=k_gw-k_bz
1778        if (Ltg_kgw%ibzq(iq_bz)==1) then
1779          bmask=.FALSE.; bmask(1:Wfd%nband(ik_ibz,spin))=.TRUE.
1780          if (PRESENT(bks_mask)) bmask = bks_mask(:,ik_bz,spin)
1781          call wfd_distribute_bands(Wfd,ik_ibz,spin,my_nband,my_band_list,got=get_more,bmask=bmask)
1782          if (my_nband>0) proc_distrb(my_band_list(1:my_nband),ik_bz,spin)=Wfd%my_rank
1783        end if
1784      end do
1785 
1786    else
1787      ! No symmetries for this spin. Divide the full BZ among procs.
1788      do ik_bz=1,Kmesh%nbz
1789        ik_ibz = Kmesh%tab(ik_bz)
1790        bmask=.FALSE.; bmask(1:Wfd%nband(ik_ibz,spin))=.TRUE.
1791        if (PRESENT(bks_mask)) bmask = bks_mask(:,ik_bz,spin)
1792        call wfd_distribute_bands(Wfd,ik_ibz,spin,my_nband,my_band_list,got=get_more,bmask=bmask)
1793        if (my_nband>0) proc_distrb(my_band_list(1:my_nband),ik_bz,spin)=Wfd%my_rank
1794      end do
1795    end if
1796  end do ! spin
1797 
1798  if (PRESENT(global)) then
1799    if (global) then ! Each node will have the same table so that it will know how the tasks are distributed.
1800      proc_distrb = proc_distrb + 1
1801      where (proc_distrb == xmpi_undefined_rank+1)
1802        proc_distrb = 0
1803      end where
1804      call xmpi_sum(proc_distrb,Wfd%comm,ierr)
1805      where (proc_distrb == 0)
1806        proc_distrb = xmpi_undefined_rank
1807      elsewhere
1808        proc_distrb = proc_distrb -1
1809      end where
1810      !where (proc_distrb /= xmpi_undefined_rank)
1811      !  ltest = (ANY(proc_distrb == (/(ii,ii=0,Wfd%nproc-1)/)))
1812      !end where
1813      !if (.not.ltest) then
1814      !  write(std_out,*)proc_distrb
1815      !  MSG_BUG("Bug in the generation of proc_distrb table")
1816      !end if
1817    end if
1818  end if
1819 
1820  my_nbks = COUNT(proc_distrb==Wfd%my_rank)
1821  if (PRESENT(got)) got=get_more
1822 
1823 end subroutine sigma_distribute_bks

m_sigma/sigma_free [ Functions ]

[ Top ] [ m_sigma ] [ Functions ]

NAME

 sigma_free

FUNCTION

  Deallocate all associated pointers defined in the sigma_t data type.

INPUTS

OUTPUT

PARENTS

      sigma

CHILDREN

      findqg0,wfd_distribute_bands,wfd_update_bkstab,xmpi_sum

SOURCE

1097 subroutine sigma_free(Sr)
1098 
1099 
1100 !This section has been created automatically by the script Abilint (TD).
1101 !Do not modify the following lines by hand.
1102 #undef ABI_FUNC
1103 #define ABI_FUNC 'sigma_free'
1104 !End of the abilint section
1105 
1106  implicit none
1107 
1108 !Arguments ------------------------------------
1109 !scalars
1110  type(sigma_t),intent(inout) :: Sr
1111 
1112 ! *************************************************************************
1113 
1114  !@sigma_t
1115 !integer
1116  if (allocated(Sr%maxbnd)) then
1117    ABI_FREE(Sr%maxbnd)
1118  end if
1119  if (allocated(Sr%minbnd)) then
1120    ABI_FREE(Sr%minbnd)
1121  end if
1122 
1123 !real
1124  if (allocated(Sr%degwgap)) then
1125    ABI_FREE(Sr%degwgap)
1126  end if
1127  if (allocated(Sr%egwgap)) then
1128    ABI_FREE(Sr%egwgap)
1129  end if
1130  if (allocated(Sr%en_qp_diago)) then
1131    ABI_FREE(Sr%en_qp_diago)
1132  end if
1133  if (allocated(Sr%e0)) then
1134    ABI_FREE(Sr%e0)
1135  end if
1136  if (allocated(Sr%e0gap)) then
1137    ABI_FREE(Sr%e0gap)
1138  end if
1139  if (allocated(Sr%omega_r)) then
1140    ABI_FREE(Sr%omega_r)
1141  end if
1142  if (allocated(Sr%kptgw)) then
1143    ABI_FREE(Sr%kptgw)
1144  end if
1145  if (allocated(Sr%sigxme)) then
1146    ABI_FREE(Sr%sigxme)
1147  end if
1148  if (allocated(Sr%x_mat)) then
1149    ABI_FREE(Sr%x_mat)
1150  end if
1151  if (allocated(Sr%vxcme)) then
1152    ABI_FREE(Sr%vxcme)
1153  end if
1154  if (allocated(Sr%vUme)) then
1155    ABI_FREE(Sr%vUme)
1156  end if
1157 
1158 !complex
1159  if (allocated(Sr%degw)) then
1160    ABI_FREE(Sr%degw)
1161  end if
1162  if (allocated(Sr%dsigmee0)) then
1163    ABI_FREE(Sr%dsigmee0)
1164  end if
1165  if (allocated(Sr%egw)) then
1166    ABI_FREE(Sr%egw)
1167  end if
1168  if (allocated(Sr%eigvec_qp)) then
1169    ABI_FREE(Sr%eigvec_qp)
1170  end if
1171  if (allocated(Sr%m_lda_to_qp)) then
1172    ABI_FREE(Sr%m_lda_to_qp)
1173  end if
1174  if (allocated(Sr%hhartree)) then
1175    ABI_FREE(Sr%hhartree)
1176  end if
1177  if (allocated(Sr%sigcme)) then
1178    ABI_FREE(Sr%sigcme)
1179  end if
1180  if (allocated(Sr%sigmee)) then
1181    ABI_FREE(Sr%sigmee)
1182  end if
1183  if (allocated(Sr%sigcmee0)) then
1184    ABI_FREE(Sr%sigcmee0)
1185  end if
1186  if (allocated(Sr%sigcmesi)) then
1187    ABI_FREE(Sr%sigcmesi)
1188  end if
1189  if (allocated(Sr%sigcme4sd)) then
1190    ABI_FREE(Sr%sigcme4sd)
1191  end if
1192  if (allocated(Sr%sigxcme)) then
1193    ABI_FREE(Sr%sigxcme)
1194  end if
1195  if (allocated(Sr%sigxcmesi)) then
1196    ABI_FREE(Sr%sigxcmesi)
1197  end if
1198  if (allocated(Sr%sigxcme4sd)) then
1199    ABI_FREE(Sr%sigxcme4sd)
1200  end if
1201  if (allocated(Sr%ze0)) then
1202    ABI_FREE(Sr%ze0)
1203  end if
1204  if (allocated(Sr%omega_i)) then
1205    ABI_FREE(Sr%omega_i)
1206  end if
1207  if (allocated(Sr%omega4sd)) then
1208    ABI_FREE(Sr%omega4sd)
1209  end if
1210 
1211 end subroutine sigma_free

m_sigma/sigma_get_exene [ Functions ]

[ Top ] [ m_sigma ] [ Functions ]

NAME

  sigma_get_exene

FUNCTION

  Compute exchange energy.

INPUTS

  sigma<sigma_t>=Sigma results
  kmesh<kmesh_t>=BZ sampling.
  bands<band_t>=Bands with occupation factors

PARENTS

SOURCE

1232 pure function sigma_get_exene(sigma,kmesh,bands) result(ex_energy)
1233 
1234 
1235 !This section has been created automatically by the script Abilint (TD).
1236 !Do not modify the following lines by hand.
1237 #undef ABI_FUNC
1238 #define ABI_FUNC 'sigma_get_exene'
1239 !End of the abilint section
1240 
1241  implicit none
1242 
1243 !Arguments ------------------------------------
1244 !scalars
1245  real(dp) :: ex_energy
1246  type(sigma_t),intent(in) :: sigma
1247  type(kmesh_t),intent(in) :: kmesh
1248  type(ebands_t),intent(in) :: bands
1249 
1250 !Local variables-------------------------------
1251 !scalars
1252  integer :: ik,ib,spin
1253  real(dp) :: wtk,occ_bks
1254 
1255 ! *************************************************************************
1256 
1257  ex_energy = zero
1258 
1259  do spin=1,sigma%nsppol
1260    do ik=1,sigma%nkibz
1261      wtk = kmesh%wt(ik)
1262      do ib=sigma%b1gw,sigma%b2gw
1263        occ_bks = bands%occ(ib,ik,spin)
1264        if (sigma%nsig_ab==1) then
1265          ex_energy = ex_energy + half * occ_bks * wtk * sigma%sigxme(ib,ik,spin)
1266        else
1267          ex_energy = ex_energy + half * occ_bks * wtk * SUM(sigma%sigxme(ib,ik,:))
1268        end if
1269      end do
1270    end do
1271  end do
1272 
1273 end function sigma_get_exene

m_sigma/sigma_init [ Functions ]

[ Top ] [ m_sigma ] [ Functions ]

NAME

 sigma_init

FUNCTION

 Main creation method for the sigma_t data type.

INPUTS

 usepawu=1 if we used LDA+U as starting point (only for PAW)

OUTPUT

TODO

  Write documentation.

PARENTS

      sigma

CHILDREN

      findqg0,wfd_distribute_bands,wfd_update_bkstab,xmpi_sum

SOURCE

 960 subroutine sigma_init(Sigp,nkibz,usepawu,Sr)
 961 
 962 
 963 !This section has been created automatically by the script Abilint (TD).
 964 !Do not modify the following lines by hand.
 965 #undef ABI_FUNC
 966 #define ABI_FUNC 'sigma_init'
 967 !End of the abilint section
 968 
 969  implicit none
 970 
 971 !Arguments ------------------------------------
 972  integer,intent(in) :: nkibz,usepawu
 973 !scalars
 974  type(sigparams_t),intent(in) :: Sigp
 975  type(sigma_t),intent(inout) :: Sr
 976 
 977 !Local variables-------------------------------
 978 !scalars
 979  integer :: b1gw,b2gw,mod10
 980 
 981 ! *************************************************************************
 982 
 983  !@sigma_t
 984  mod10=MOD(Sigp%gwcalctyp,10)
 985 
 986  ! Copy important dimensions
 987  Sr%nkptgw     =Sigp%nkptgw
 988  Sr%gwcalctyp  =Sigp%gwcalctyp
 989  Sr%deltae     =Sigp%deltae
 990  Sr%maxomega4sd=Sigp%maxomega4sd
 991  Sr%maxomega_r =Sigp%maxomega_r
 992  Sr%scissor_ene=Sigp%mbpt_sciss
 993 
 994  !FIXME this should be done in sigma_allocate
 995  ABI_MALLOC(Sr%minbnd,(Sr%nkptgw,Sigp%nsppol))
 996  ABI_MALLOC(Sr%maxbnd,(Sr%nkptgw,Sigp%nsppol))
 997  Sr%minbnd=Sigp%minbnd; Sr%maxbnd=Sigp%maxbnd
 998  ABI_MALLOC(Sr%kptgw,(3,Sr%nkptgw))
 999  Sr%kptgw=Sigp%kptgw
1000 
1001  Sr%b1gw     =Sigp%minbdgw ! * min and Max GW band index over k and spin.
1002  Sr%b2gw     =Sigp%maxbdgw !   Used to dimension arrays.
1003  Sr%nbnds    =Sigp%nbnds
1004  Sr%nkibz    =nkibz
1005  Sr%nsppol   =Sigp%nsppol
1006  Sr%nsig_ab  =Sigp%nsig_ab
1007  Sr%nomega_r =Sigp%nomegasr  !FIXME change name
1008  Sr%nomega_i =Sigp%nomegasi
1009  Sr%nomega4sd=Sigp%nomegasrd
1010  Sr%usepawu  =usepawu
1011 
1012  !======================================================
1013  ! === Allocate arrays in the sigma_t datatype ===
1014  !======================================================
1015  b1gw=Sr%b1gw
1016  b2gw=Sr%b2gw
1017 
1018  !TODO write routine to allocate all this stuff
1019 
1020  ! hhartree(b1,b2,k,s)= <b1,k,s|T+v_{loc}+v_{nl}+v_{H}|b2,k,s>
1021  ABI_CALLOC(Sr%hhartree,(b1gw:b2gw,b1gw:b2gw,Sr%nkibz,Sr%nsppol*Sr%nsig_ab))
1022 
1023  ! QP amplitudes and energies
1024  ABI_CALLOC(Sr%en_qp_diago,(Sr%nbnds,Sr%nkibz,Sr%nsppol))
1025  ABI_CALLOC(Sr%eigvec_qp,(Sr%nbnds,Sr%nbnds,Sr%nkibz,Sr%nsppol))
1026 
1027  ! Dont know if it is better to do this here or in the sigma
1028  ! * Initialize with KS wavefunctions and energies
1029  !do ib=1,Sr%nbnds
1030  ! Sr%en_qp_diago(ib,:,:)=en(:,ib,:)
1031  ! Sr%eigvec_qp(ib,ib,:,:)=cone
1032  !end do
1033 
1034  ABI_CALLOC(Sr%vxcme,(b1gw:b2gw,Sr%nkibz,Sr%nsppol*Sr%nsig_ab))
1035  ABI_CALLOC(Sr%vUme,(b1gw:b2gw,Sr%nkibz,Sr%nsppol*Sr%nsig_ab))
1036  ABI_CALLOC(Sr%sigxme,(b1gw:b2gw,Sr%nkibz,Sr%nsppol*Sr%nsig_ab))
1037  ABI_CALLOC(Sr%x_mat,(b1gw:b2gw,b1gw:b2gw,Sr%nkibz,Sr%nsppol*Sr%nsig_ab))
1038  ABI_CALLOC(Sr%sigcme,(b1gw:b2gw,Sr%nkibz,Sr%nomega_r,Sr%nsppol*Sr%nsig_ab))
1039  ABI_CALLOC(Sr%sigxcme,(b1gw:b2gw,Sr%nkibz,Sr%nomega_r,Sr%nsppol*Sr%nsig_ab))
1040  ABI_CALLOC(Sr%sigcmee0,(b1gw:b2gw,Sr%nkibz,Sr%nsppol*Sr%nsig_ab))
1041  ABI_CALLOC(Sr%ze0,(b1gw:b2gw,Sr%nkibz,Sr%nsppol))
1042  ABI_CALLOC(Sr%dsigmee0,(b1gw:b2gw,Sr%nkibz,Sr%nsppol*Sr%nsig_ab))
1043  ABI_CALLOC(Sr%sigmee,(b1gw:b2gw,Sr%nkibz,Sr%nsppol*Sr%nsig_ab))
1044  ABI_CALLOC(Sr%degw,(b1gw:b2gw,Sr%nkibz,Sr%nsppol))
1045  ABI_CALLOC(Sr%e0,(Sr%nbnds,Sr%nkibz,Sr%nsppol))
1046  ABI_CALLOC(Sr%egw,(Sr%nbnds,Sr%nkibz,Sr%nsppol))
1047  ABI_CALLOC(Sr%e0gap,(Sr%nkibz,Sr%nsppol))
1048  ABI_CALLOC(Sr%degwgap,(Sr%nkibz,Sr%nsppol))
1049  ABI_CALLOC(Sr%egwgap,(Sr%nkibz,Sr%nsppol))
1050 
1051  ! These quantities are used to evaluate $\Sigma(E)$ around the KS\QP eigenvalue
1052  ABI_CALLOC(Sr%omega4sd,(b1gw:b2gw,Sr%nkibz,Sr%nomega4sd,Sr%nsppol))
1053  ABI_CALLOC(Sr%sigcme4sd,(b1gw:b2gw,Sr%nkibz,Sr%nomega4sd,Sr%nsppol*Sr%nsig_ab))
1054  ABI_CALLOC(Sr%sigxcme4sd,(b1gw:b2gw,Sr%nkibz,Sr%nomega4sd,Sr%nsppol*Sr%nsig_ab))
1055 
1056  !TODO Find  better treatment
1057  ! Mesh along the real axis.
1058  if (Sr%nomega_r>0) then
1059    ABI_MALLOC(Sr%omega_r,(Sr%nomega_r))
1060    Sr%omega_r(:)=Sigp%omega_r(:)
1061  end if
1062 
1063  ! Analytical Continuation
1064  ! FIXME omegasi should not be in Sigp% here we should construct the mesh
1065  if (mod10==1) then
1066    ABI_MALLOC(Sr%omega_i,(Sr%nomega_i))
1067    Sr%omega_i=Sigp%omegasi
1068    ABI_MALLOC(Sr%sigcmesi ,(b1gw:b2gw,Sr%nkibz,Sr%nomega_i,Sr%nsppol*Sr%nsig_ab))
1069    ABI_MALLOC(Sr%sigxcmesi,(b1gw:b2gw,Sr%nkibz,Sr%nomega_i,Sr%nsppol*Sr%nsig_ab))
1070    Sr%sigcmesi=czero; Sr%sigxcmesi=czero
1071  end if
1072 
1073 end subroutine sigma_init

m_sigma/sigma_ncwrite [ Functions ]

[ Top ] [ m_sigma ] [ Functions ]

NAME

 sigma_ncwrite

FUNCTION

  Save the data stored in the sigma_t data type on a NETCDF file.

INPUTS

  filename

OUTPUT

PARENTS

      sigma

CHILDREN

SOURCE

1418 integer function sigma_ncwrite(Sigp,Er,Sr,ncid) result (ncerr)
1419 
1420 
1421 !This section has been created automatically by the script Abilint (TD).
1422 !Do not modify the following lines by hand.
1423 #undef ABI_FUNC
1424 #define ABI_FUNC 'sigma_ncwrite'
1425 !End of the abilint section
1426 
1427  implicit none
1428 
1429 !Arguments ------------------------------------
1430 !scalars
1431  integer,intent(in) :: ncid
1432  type(sigparams_t),target,intent(in) :: Sigp
1433  type(Epsilonm1_results),target,intent(in) :: Er
1434  type(sigma_t),target,intent(in) :: Sr
1435 
1436 !Local variables ---------------------------------------
1437 #ifdef HAVE_NETCDF
1438 !scalars
1439  integer :: nbgw,ndim_sig,b1gw,b2gw,cplex
1440  !character(len=500) :: msg
1441 !arrays
1442  real(dp),allocatable :: rdata2(:,:),rdata4(:,:,:,:),rdata5(:,:,:,:,:)
1443 
1444 ! *************************************************************************
1445 
1446  !@sigma_t
1447  cplex=2; b1gw=Sr%b1gw; b2gw=Sr%b2gw; nbgw=b2gw-b1gw+1
1448  ndim_sig=Sr%nsppol*Sr%nsig_ab
1449 
1450  ncerr = nctk_def_dims(ncid, [nctkdim_t("cplex", cplex), nctkdim_t("b1gw", sr%b1gw), nctkdim_t("b2gw", sr%b2gw),&
1451    nctkdim_t("nbgw", nbgw), nctkdim_t("nkptgw", sr%nkptgw), nctkdim_t("ndim_sig", ndim_sig), &
1452    nctkdim_t("nomega4sd", sr%nomega4sd), nctkdim_t("nsig_ab", sr%nsig_ab)], defmode=.True.)
1453  NCF_CHECK(ncerr)
1454 
1455  ! No. of real frequencies, might be zero.
1456  if (Sr%nomega_r>0) then
1457    NCF_CHECK(nctk_def_dims(ncid, nctkdim_t("nomega_r", Sr%nomega_r)))
1458  end if
1459 
1460  ! No. of imaginary frequencies, might be zero.
1461  if (Sr%nomega_i>0) then
1462    NCF_CHECK(nctk_def_dims(ncid, nctkdim_t("nomega_i", Sr%nomega_i)))
1463  end if
1464 
1465  ! =======================
1466  ! == Define variables ===
1467  ! =======================
1468  ! parameters of the calculation.
1469 
1470  ncerr = nctk_def_iscalars(ncid, [character(len=nctk_slen) :: 'sigma_nband', 'scr_nband', 'gwcalctyp', 'usepawu'])
1471  NCF_CHECK(ncerr)
1472 
1473  ncerr = nctk_def_dpscalars(ncid, [character(len=nctk_slen) :: &
1474 &  'ecutwfn', 'ecuteps', 'ecutsigx', 'omegasrdmax', 'deltae', 'omegasrmax', 'scissor_ene'])
1475  NCF_CHECK(ncerr)
1476 
1477  ! TODO: Decrease size of file: Remove arrays whose size scale as mband ** 2
1478  ! especially those that are not commonly used e.g. hhartree.
1479  ncerr = nctk_def_arrays(ncid, [&
1480    nctkarr_t("kptgw", "dp", "number_of_reduced_dimensions, nkptgw"),&
1481    nctkarr_t("minbnd", "i", "nkptgw, number_of_spins"),&
1482    nctkarr_t("maxbnd", "i", "nkptgw, number_of_spins"), &
1483    nctkarr_t('degwgap', "dp", 'number_of_kpoints, number_of_spins'),&
1484    nctkarr_t('egwgap', "dp", 'number_of_kpoints, number_of_spins'),&
1485    nctkarr_t('en_qp_diago', "dp",'max_number_of_states, number_of_kpoints, number_of_spins'),&
1486    nctkarr_t('e0', "dp", 'max_number_of_states, number_of_kpoints, number_of_spins'),&
1487    nctkarr_t('e0gap', "dp", 'number_of_kpoints, number_of_spins'),&
1488    nctkarr_t('sigxme', "dp", 'nbgw, number_of_kpoints, ndim_sig'),&
1489    nctkarr_t('vxcme', "dp", 'nbgw, number_of_kpoints, ndim_sig'),&
1490    nctkarr_t('degw', "dp", 'cplex, nbgw, number_of_kpoints, number_of_spins'),&
1491    nctkarr_t('dsigmee0', "dp", 'cplex, nbgw, number_of_kpoints, ndim_sig'),&
1492    nctkarr_t('egw', "dp",'cplex, max_number_of_states, number_of_kpoints, number_of_spins'),&
1493    nctkarr_t('eigvec_qp', "dp",'cplex, max_number_of_states, max_number_of_states, number_of_kpoints, number_of_spins'),&
1494    nctkarr_t('hhartree', "dp",'cplex, nbgw, nbgw, number_of_kpoints, ndim_sig'),&
1495    nctkarr_t('sigmee', "dp", 'cplex, nbgw, number_of_kpoints, ndim_sig'),&
1496    nctkarr_t('sigcmee0', "dp",'cplex, nbgw, number_of_kpoints, ndim_sig'),&
1497    nctkarr_t('sigcmesi', "dp",'cplex, nbgw, number_of_kpoints, ndim_sig'),&
1498    nctkarr_t('sigcme4sd', "dp",'cplex, nbgw, number_of_kpoints, nomega4sd, ndim_sig'),&
1499    nctkarr_t('sigxcme4sd', "dp", 'cplex, nbgw, number_of_kpoints, nomega4sd, ndim_sig'),&
1500    nctkarr_t('ze0',"dp", 'cplex, nbgw, number_of_kpoints, number_of_spins'),&
1501    nctkarr_t('omega4sd', "dp", 'cplex, nbgw, number_of_kpoints, nomega4sd, number_of_spins')])
1502  NCF_CHECK(ncerr)
1503 
1504  if (Sr%usepawu==0) then
1505    ncerr = nctk_def_arrays(ncid, nctkarr_t("vUme", "dp", 'nbgw, number_of_kpoints, ndim_sig'))
1506    NCF_CHECK(ncerr)
1507  end if
1508 
1509  if (Sr%nomega_r>0) then
1510    ncerr = nctk_def_arrays(ncid, [&
1511      nctkarr_t('omega_r', "dp", "nomega_r"),&
1512      nctkarr_t('sigcme', "dp", 'cplex, nbgw, number_of_kpoints, nomega_r, ndim_sig'),&
1513      nctkarr_t('sigxcme', "dp", 'cplex, nbgw, number_of_kpoints, nomega_r, ndim_sig')])
1514    NCF_CHECK(ncerr)
1515  end if
1516 
1517  if (Sr%nomega_i>0) then
1518    ncerr = nctk_def_arrays(ncid, [&
1519      nctkarr_t('sigxcmesi', "dp", 'cplex, nbgw, number_of_kpoints, nomega_i, ndim_sig'),&
1520      nctkarr_t('omega_i', "dp", 'cplex, nomega_i')])
1521    NCF_CHECK(ncerr)
1522  end if
1523 
1524  if (allocated(sr%m_lda_to_qp)) then
1525    ncerr = nctk_def_arrays(ncid, [nctkarr_t('m_lda_to_qp', "dp", &
1526 &       "cplex, max_number_of_states, max_number_of_states, number_of_kpoints, number_of_spins")])
1527    NCF_CHECK(ncerr)
1528  end if
1529 
1530  ! =====================
1531  ! === Start writing ===
1532  ! =====================
1533  ! vid is a small helper function to get the varid from ncid
1534  NCF_CHECK(nctk_set_datamode(ncid))
1535  NCF_CHECK(nf90_put_var(ncid, vid('ecutwfn'), Sigp%ecutwfn))
1536  NCF_CHECK(nf90_put_var(ncid, vid('ecuteps'), Sigp%ecuteps))
1537  NCF_CHECK(nf90_put_var(ncid, vid('ecutsigx'), Sigp%ecutsigx))
1538  NCF_CHECK(nf90_put_var(ncid, vid('sigma_nband'), Sigp%nbnds))
1539  NCF_CHECK(nf90_put_var(ncid, vid('scr_nband'), Er%Hscr%nbnds_used))
1540  NCF_CHECK(nf90_put_var(ncid, vid('gwcalctyp'), Sr%gwcalctyp))
1541  NCF_CHECK(nf90_put_var(ncid, vid('usepawu'), Sr%usepawu))
1542  NCF_CHECK(nf90_put_var(ncid, vid('kptgw'), Sr%kptgw))
1543  NCF_CHECK(nf90_put_var(ncid, vid('minbnd'), Sr%minbnd))
1544  NCF_CHECK(nf90_put_var(ncid, vid('maxbnd'),Sr%maxbnd))
1545 
1546  NCF_CHECK(nf90_put_var(ncid, vid('omegasrdmax'), Sr%maxomega4sd*Ha_eV))
1547  NCF_CHECK(nf90_put_var(ncid, vid('deltae'), Sr%deltae*Ha_eV))
1548  NCF_CHECK(nf90_put_var(ncid, vid('omegasrmax'), Sr%maxomega_r*Ha_eV))
1549  NCF_CHECK(nf90_put_var(ncid, vid('scissor_ene'), Sr%scissor_ene*Ha_eV))
1550  NCF_CHECK(nf90_put_var(ncid, vid('degwgap'), Sr%degwgap*Ha_eV))
1551  NCF_CHECK(nf90_put_var(ncid, vid('egwgap'), Sr%egwgap*Ha_eV))
1552  NCF_CHECK(nf90_put_var(ncid, vid('en_qp_diago'), Sr%en_qp_diago*Ha_eV))
1553  NCF_CHECK(nf90_put_var(ncid, vid('e0'), Sr%e0*Ha_eV))
1554  NCF_CHECK(nf90_put_var(ncid, vid('e0gap'), Sr%e0gap*Ha_eV))
1555 
1556  if (Sr%nomega_r>0) then
1557    NCF_CHECK(nf90_put_var(ncid, vid('omega_r'), Sr%omega_r*Ha_eV))
1558  end if
1559 
1560  NCF_CHECK(nf90_put_var(ncid, vid('sigxme'), Sr%sigxme*Ha_eV))
1561  NCF_CHECK(nf90_put_var(ncid, vid('vxcme'), Sr%vxcme*Ha_eV))
1562  NCF_CHECK(nf90_put_var(ncid, vid('vUme'), Sr%vUme*Ha_eV))
1563 
1564  ! * Have to transfer complex arrays
1565  ABI_MALLOC(rdata4,(cplex,b1gw:b2gw,Sr%nkibz,Sr%nsppol))
1566  rdata4=c2r(Sr%degw)
1567  NCF_CHECK(nf90_put_var(ncid, vid('degw'), rdata4*Ha_eV))
1568  ABI_FREE(rdata4)
1569 
1570  ABI_MALLOC(rdata4,(cplex,b1gw:b2gw,Sr%nkibz,Sr%nsppol*Sr%nsig_ab))
1571  rdata4=c2r(Sr%dsigmee0)
1572  NCF_CHECK(nf90_put_var(ncid, vid('dsigmee0'), rdata4))
1573  ABI_FREE(rdata4)
1574 
1575  ABI_MALLOC(rdata4,(cplex,Sr%nbnds,Sr%nkibz,Sr%nsppol))
1576  rdata4=c2r(Sr%egw)
1577  NCF_CHECK(nf90_put_var(ncid, vid('egw'), rdata4*Ha_eV))
1578  ABI_FREE(rdata4)
1579 
1580  ABI_MALLOC(rdata5,(cplex,Sr%nbnds,Sr%nbnds,Sr%nkibz,Sr%nsppol))
1581  rdata5=c2r(Sr%eigvec_qp)
1582  NCF_CHECK(nf90_put_var(ncid, vid('eigvec_qp'), rdata5))
1583  ABI_FREE(rdata5)
1584 
1585  ABI_MALLOC(rdata5,(cplex,nbgw,nbgw,Sr%nkibz,Sr%nsppol*Sr%nsig_ab))
1586  rdata5=c2r(Sr%hhartree)
1587  NCF_CHECK(nf90_put_var(ncid, vid('hhartree'), rdata5*Ha_eV))
1588  ABI_FREE(rdata5)
1589 
1590  if (Sr%nomega_r>0) then
1591    ABI_MALLOC(rdata5,(cplex,nbgw,Sr%nkibz,Sr%nomega_r,Sr%nsppol*Sr%nsig_ab))
1592    rdata5=c2r(Sr%sigcme)
1593    NCF_CHECK(nf90_put_var(ncid, vid('sigcme'), rdata5*Ha_eV))
1594    ABI_FREE(rdata5)
1595  end if
1596 
1597  ABI_MALLOC(rdata4,(cplex,nbgw,Sr%nkibz,Sr%nsppol*Sr%nsig_ab))
1598  rdata4=c2r(Sr%sigmee)
1599  NCF_CHECK(nf90_put_var(ncid, vid('sigmee'), rdata4*Ha_eV))
1600  ABI_FREE(rdata4)
1601 
1602  ABI_MALLOC(rdata4,(cplex,nbgw,Sr%nkibz,Sr%nsppol*Sr%nsig_ab))
1603  rdata4=c2r(Sr%sigcmee0)
1604  NCF_CHECK(nf90_put_var(ncid, vid('sigcmee0'), rdata4*Ha_eV))
1605  ABI_FREE(rdata4)
1606 
1607  if (Sr%nomega_i>0) then
1608   ABI_MALLOC(rdata5,(cplex,nbgw,Sr%nkibz,Sr%nomega_i,Sr%nsppol*Sr%nsig_ab))
1609   rdata5=c2r(Sr%sigcmesi)
1610   NCF_CHECK(nf90_put_var(ncid, vid('sigcmesi'), rdata5*Ha_eV))
1611   ABI_FREE(rdata5)
1612  end if
1613 
1614  ABI_MALLOC(rdata5,(cplex,nbgw,Sr%nkibz,Sr%nomega4sd,Sr%nsppol*Sr%nsig_ab))
1615  rdata5=c2r(Sr%sigcme4sd)
1616  NCF_CHECK(nf90_put_var(ncid, vid('sigcme4sd'), rdata5*Ha_eV))
1617  ABI_FREE(rdata5)
1618 
1619  if (Sr%nomega_r>0) then
1620    ABI_MALLOC(rdata5,(cplex,nbgw,Sr%nkibz,Sr%nomega_r,Sr%nsppol*Sr%nsig_ab))
1621    rdata5=c2r(Sr%sigxcme)
1622    NCF_CHECK(nf90_put_var(ncid, vid('sigxcme'), rdata5*Ha_eV))
1623    ABI_FREE(rdata5)
1624  end if
1625 
1626  if (Sr%nomega_i>0) then
1627    ABI_MALLOC(rdata5,(cplex,nbgw,Sr%nkibz,Sr%nomega_i,Sr%nsppol*Sr%nsig_ab))
1628    rdata5=c2r(Sr%sigxcmesi)
1629    NCF_CHECK(nf90_put_var(ncid, vid('sigxcmesi'), rdata5*Ha_eV))
1630    ABI_FREE(rdata5)
1631  end if
1632 
1633  if (allocated(sr%m_lda_to_qp)) then
1634    ABI_MALLOC(rdata5,(cplex,Sr%nbnds,Sr%nbnds,Sr%nkibz,Sr%nsppol))
1635    rdata5=c2r(Sr%m_lda_to_qp)
1636    NCF_CHECK(nf90_put_var(ncid, vid('m_lda_to_qp'), rdata5))
1637    ABI_FREE(rdata5)
1638  end if
1639 
1640  ABI_MALLOC(rdata5,(cplex,nbgw,Sr%nkibz,Sr%nomega4sd,Sr%nsppol*Sr%nsig_ab))
1641  rdata5=c2r(Sr%sigxcme4sd)
1642  NCF_CHECK(nf90_put_var(ncid, vid('sigxcme4sd'), rdata5*Ha_eV))
1643  ABI_FREE(rdata5)
1644 
1645  ABI_MALLOC(rdata4,(cplex,nbgw,Sr%nkibz,Sr%nsppol))
1646  rdata4=c2r(Sr%ze0)
1647  NCF_CHECK(nf90_put_var(ncid, vid('ze0'), rdata4))
1648  ABI_FREE(rdata4)
1649 
1650  if (Sr%nomega_i>0) then
1651    ABI_MALLOC(rdata2,(cplex,Sr%nomega_i))
1652    rdata2=c2r(Sr%omega_i)
1653    NCF_CHECK(nf90_put_var(ncid, vid('omega_i'), rdata2*Ha_eV))
1654    ABI_FREE(rdata2)
1655  end if
1656 
1657  ABI_MALLOC(rdata5,(cplex,nbgw,Sr%nkibz,Sr%nomega4sd,Sr%nsppol))
1658  rdata5=c2r(Sr%omega4sd)
1659  NCF_CHECK(nf90_put_var(ncid, vid('omega4sd'), rdata5*Ha_eV))
1660  ABI_FREE(rdata5)
1661 
1662 #else
1663   MSG_ERROR('netcdf support is not activated.')
1664 #endif
1665 
1666 contains
1667  integer function vid(vname)
1668 
1669 
1670 !This section has been created automatically by the script Abilint (TD).
1671 !Do not modify the following lines by hand.
1672 #undef ABI_FUNC
1673 #define ABI_FUNC 'vid'
1674 !End of the abilint section
1675 
1676    character(len=*),intent(in) :: vname
1677    vid = nctk_idname(ncid, vname)
1678  end function vid
1679 
1680 end function sigma_ncwrite

m_sigma/sigma_t [ Types ]

[ Top ] [ m_sigma ] [ Types ]

NAME

 sigma_t

FUNCTION

 For the GW part of ABINIT, the sigma_t structured datatype
 gather the results of a sigma calculation.

TODO

   ragged arrays (nk,nsppol) --> values ?

SOURCE

 68  type,public :: sigma_t
 69 
 70   integer :: b1gw,b2gw      ! min and Max gw band indeces over spin and k-points (used to dimension arrays)
 71   integer :: gwcalctyp      ! Flag defining the calculation type.
 72   integer :: nkptgw         ! No. of points calculated
 73   integer :: nkibz          ! No. of irreducible k-points.
 74   integer :: nbnds          ! Total number of bands
 75   integer :: nomega_r       ! No. of real frequencies for the spectral function.
 76   integer :: nomega_i       ! No. of frequencies along the imaginary axis.
 77   integer :: nomega4sd      ! No. of real frequencies to evaluate the derivative of $\Sigma(E)$.
 78   integer :: nsig_ab        ! 1 if nspinor=1,4 for noncollinear case.
 79   integer :: nsppol         ! No. of spin polarizations.
 80   integer :: usepawu        ! 1 if we are using LDA+U as starting point (only for PAW)
 81 
 82   real(dp) :: deltae       ! Frequency step for the calculation of d\Sigma/dE
 83   real(dp) :: maxomega4sd  ! Max frequency around E_ks for d\Sigma/dE.
 84   real(dp) :: maxomega_r   ! Max frequency for spectral function.
 85   real(dp) :: scissor_ene  ! Scissor energy value. zero for None.
 86 
 87   integer,allocatable :: maxbnd(:,:)
 88   ! maxbnd(nkptgw,nsppol)
 89   ! Max band index considered in GW for this k-point.
 90 
 91   integer,allocatable :: minbnd(:,:)
 92   ! minbnd(nkptgw,nsppol)
 93   ! Min band index considered in GW for this k-point.
 94 
 95   !real(dp),allocatable :: ame(:,:,:)
 96   ! ame(nbnds,nkibz,nomega))
 97   ! Diagonal matrix elements of the spectral function.
 98   ! Commented out, it can be calculated from the other quantities
 99 
100   real(dp),allocatable :: degwgap(:,:)
101   ! degwgap(nkibz,nsppol)
102   ! Difference btw the QP and the KS optical gap.
103 
104   real(dp),allocatable :: egwgap(:,:)
105   ! egwgap(nkibz,nsppol))
106   ! QP optical gap at each k-point and spin.
107 
108   real(dp),allocatable :: en_qp_diago(:,:,:)
109   ! en_qp_diago(nbnds,nkibz,nsppol))
110   ! QP energies obtained from the diagonalization of the Hermitian approximation to Sigma (QPSCGW)
111 
112   real(dp),allocatable :: e0(:,:,:)
113   ! e0(nbnds,nkibz,nsppol)
114   ! KS eigenvalues for each band, k-point and spin. In case of self-consistent?
115 
116   real(dp),allocatable :: e0gap(:,:)
117   ! e0gap(nkibz,nsppol),
118   ! KS gap at each k-point, for each spin.
119 
120   real(dp),allocatable :: omega_r(:)
121   ! omega_r(nomega_r)
122   ! real frequencies used for the self energy.
123 
124   real(dp),allocatable :: kptgw(:,:)
125   ! kptgw(3,nkptgw)
126   ! ! TODO there is a similar array in sigparams_t
127   ! List of calculated k-points.
128 
129   real(dp),allocatable :: sigxme(:,:,:)
130   ! sigxme(b1gw:b2gw,nkibz,nsppol*nsig_ab))
131   ! Diagonal matrix elements $\<nks|\Sigma_x|nks\>$
132 
133   complex(dp),allocatable :: x_mat(:,:,:,:)
134   ! x_mat(b1gw:b2gw,b1gw:b2gw,nkibz,nsppol*nsig_ab))
135   ! Matrix elements of $\<nks|\Sigma_x|nk's\>$
136 
137   real(dp),allocatable :: vxcme(:,:,:)
138   ! vxcme(b1gw:b2gw,nkibz,nsppol*nsig_ab))
139   ! $\<nks|v_{xc}[n_val]|nks\>$ matrix elements of vxc (valence-only contribution).
140 
141   real(dp),allocatable :: vUme(:,:,:)
142   ! vUme(b1gw:b2gw,nkibz,nsppol*nsig_ab))
143   ! $\<nks|v_{U}|nks\>$ for LDA+U.
144 
145   complex(dpc),allocatable :: degw(:,:,:)
146   ! degw(b1gw:b2gw,nkibz,nsppol))
147   ! Difference between the QP and the KS energies.
148 
149   complex(dpc),allocatable :: dsigmee0(:,:,:)
150   ! dsigmee0(b1gw:b2gw,nkibz,nsppol*nsig_ab))
151   ! Derivative of $\Sigma_c(E)$ calculated at the KS eigenvalue.
152 
153   complex(dpc),allocatable :: egw(:,:,:)
154   ! degw(nbnds,nkibz,nsppol))
155   ! QP energies, $\epsilon_{nks}^{QP}$.
156 
157   complex(dpc),allocatable :: eigvec_qp(:,:,:,:)
158   ! eigvec_qp(nbnds,nbnds,nkibz,nsppol))
159   ! Expansion of the QP amplitudes in the QP basis set of the previous iteration.
160 
161   complex(dpc),allocatable :: m_lda_to_qp(:,:,:,:)
162   ! (%nbnds,%nbnds,%nibz,%nsppol))
163   !  m_lda_to_qp(ib,jb,k,s) := <\psi_{ib,k,s}^{KS}|\psi_{jb,k,s}^{QP}>
164 
165   complex(dpc),allocatable :: hhartree(:,:,:,:)
166   ! hhartree(b1gw:b2gw,b1gw:b2gw,nkibz,nsppol*nsig_ab)
167   ! $\<nks|T+v_H+v_{loc}+v_{nl}|mks\>$
168 
169   complex(dpc),allocatable :: sigcme(:,:,:,:)
170   ! sigcme(b1gw:b2gw,nkibz,nomega_r,nsppol*nsig_ab))
171   ! $\<nks|\Sigma_{c}(E)|nks\>$ at each nomega_r frequency
172 
173   complex(dpc),allocatable :: sigmee(:,:,:)
174   ! sigmee(b1gw:b2gw,nkibz,nsppol*nsig_ab))
175   ! $\Sigma_{xc}E_{KS} + (E_{QP}- E_{KS})*dSigma/dE_KS
176 
177   complex(dpc),allocatable :: sigcmee0(:,:,:)
178   ! sigcmee0(b1gw:b2gw,nkibz,nsppol*nsig_ab))
179   ! Diagonal mat. elements of $\Sigma_c(E)$ calculated at the KS energy $E_{KS}$
180 
181   complex(dpc),allocatable :: sigcmesi(:,:,:,:)
182   ! sigcmesi(b1gw:b2gw,nkibz,nomega_i,nsppol*nsig_ab))
183   ! Matrix elements of $\Sigma_c$ along the imaginary axis.
184   ! Only used in case of analytical continuation.
185 
186   complex(dpc),allocatable :: sigcme4sd(:,:,:,:)
187   ! sigcme4sd(b1gw:b2gw,nkibz,nomega4sd,nsppol*nsig_ab))
188   ! Diagonal matrix elements of \Sigma_c around the zeroth order eigenvalue (usually KS).
189 
190   complex(dpc),allocatable :: sigxcme(:,:,:,:)
191   ! sigxcme(b1gw:b2gw,nkibz,nomega_r,nsppol*nsig_ab))
192   ! $\<nks|\Sigma_{xc}(E)|nks\>$ at each real frequency frequency.
193 
194   complex(dpc),allocatable :: sigxcmesi(:,:,:,:)
195   ! sigxcmesi(b1gw:b2gw,nkibz,nomega_i,nsppol*nsig_ab))
196   ! Matrix elements of $\Sigma_{xc}$ along the imaginary axis.
197   ! Only used in case of analytical continuation.
198 
199   complex(dpc),allocatable :: sigxcme4sd(:,:,:,:)
200   ! sigxcme4sd(b1gw:b2gw,nkibz,nomega4sd,nsppol*nsig_ab))
201   ! Diagonal matrix elements of \Sigma_xc for frequencies around the zeroth order eigenvalues.
202 
203   complex(dpc),allocatable :: ze0(:,:,:)
204   ! ze0(b1gw:b2gw,nkibz,nsppol))
205   ! renormalization factor. $(1-\dfrac{\partial\Sigma_c} {\partial E_{KS}})^{-1}$
206 
207   complex(dpc),allocatable :: omega_i(:)
208   ! omegasi(nomega_i)
209   ! Frequencies along the imaginary axis used for the analytical continuation.
210 
211   complex(dpc),allocatable :: omega4sd(:,:,:,:)
212   ! omega4sd(b1gw:b2gw,nkibz,nomega4sd,nsppol).
213   ! Frequencies used to evaluate the Derivative of Sigma.
214 
215  end type sigma_t
216 
217  public  :: sigma_init                  ! Initialize the object
218  public  :: sigma_free                  ! Deallocate memory
219  public  :: sigma_get_exene             ! Compute exchange energy.
220  public  :: sigma_ncwrite               ! Write data in netcdf format.
221  public  :: write_sigma_header
222  public  :: write_sigma_results
223  public  :: print_Sigma_perturbative
224  public  :: print_Sigma_QPSC
225  public  :: sigma_distribute_bks

m_sigma/write_sigma_header [ Functions ]

[ Top ] [ m_sigma ] [ Functions ]

NAME

 write_sigma_header

FUNCTION

  Write basic info and dimensions used during the calculation
  of the QP correctoions (optdriver==4).

INPUTS

  Sigp=sigparams_t
  Cryst<crystal_t>= Info on the Crystal structure
  Kmesh<kmesh_t>= Description of the BZ sampling.

OUTPUT

  (for writing routines, no output) otherwise, should be described

NOTES

PARENTS

      sigma

CHILDREN

      findqg0,wfd_distribute_bands,wfd_update_bkstab,xmpi_sum

SOURCE

258 subroutine write_sigma_header(Sigp,Er,Cryst,Kmesh,Qmesh)
259 
260 
261 !This section has been created automatically by the script Abilint (TD).
262 !Do not modify the following lines by hand.
263 #undef ABI_FUNC
264 #define ABI_FUNC 'write_sigma_header'
265 !End of the abilint section
266 
267  implicit none
268 
269 !Arguments ------------------------------------
270 !scalars
271  type(kmesh_t),intent(in) :: Kmesh,Qmesh
272  type(crystal_t),intent(in) :: Cryst
273  type(Epsilonm1_results),intent(in) :: Er
274  type(sigparams_t),intent(in) :: Sigp
275 
276 !Local variables-------------------------------
277 !scalars
278  integer :: gwcalctyp,mod10
279  character(len=500) :: msg
280 
281 ! *************************************************************************
282 
283  msg = ' SIGMA fundamental parameters:'
284  call wrtout(std_out,msg,'COLL')
285  call wrtout(ab_out,msg,'COLL')
286 
287  gwcalctyp=Sigp%gwcalctyp
288  mod10=MOD(Sigp%gwcalctyp,10)
289 
290  SELECT CASE (mod10)
291  CASE (0)
292    write(msg,'(a,i2)')' PLASMON POLE MODEL ',Sigp%ppmodel
293  CASE (1)
294    write(msg,'(a)')' ANALYTIC CONTINUATION'
295  CASE (2)
296    write(msg,'(a)')' CONTOUR DEFORMATION'
297  CASE (5)
298    write(msg,'(a)')' Hartree-Fock'
299  CASE (6)
300    write(msg,'(a)')' Screened Exchange'
301  CASE (7)
302    write(msg,'(a)')' COHSEX'
303  CASE (8)
304    write(msg,'(a,i2)')' MODEL GW with PLASMON POLE MODEL ',Sigp%ppmodel
305  CASE (9)
306    write(msg,'(a)')' MODEL GW without PLASMON POLE MODEL'
307  CASE DEFAULT
308    write(msg,'(a,i3)')' Wrong value for Sigp%gwcalctyp = ',Sigp%gwcalctyp
309    MSG_BUG(msg)
310  END SELECT
311  call wrtout(std_out,msg,'COLL')
312  call wrtout(ab_out,msg,'COLL')
313 
314  write(msg,'(a,i12)')' number of plane-waves for SigmaX         ',Sigp%npwx
315  call wrtout(std_out,msg,'COLL')
316  call wrtout(ab_out,msg,'COLL')
317  write(msg,'(a,i12)')' number of plane-waves for SigmaC and W   ',Sigp%npwc
318  call wrtout(std_out,msg,'COLL')
319  call wrtout(ab_out,msg,'COLL')
320  write(msg,'(a,i12)')' number of plane-waves for wavefunctions  ',Sigp%npwwfn
321  call wrtout(std_out,msg,'COLL')
322  call wrtout(ab_out,msg,'COLL')
323  write(msg,'(a,i12)')' number of bands                          ',Sigp%nbnds
324  call wrtout(std_out,msg,'COLL')
325  call wrtout(ab_out,msg,'COLL')
326  write(msg,'(a,i12)')' number of independent spin polarizations ',Sigp%nsppol
327  call wrtout(std_out,msg,'COLL')
328  call wrtout(ab_out,msg,'COLL')
329  write(msg,'(a,i12)')' number of spinorial components           ',Sigp%nspinor
330  call wrtout(std_out,msg,'COLL')
331  call wrtout(ab_out,msg,'COLL')
332  write(msg,'(a,i12)')' number of k-points in IBZ                ',Kmesh%nibz
333  call wrtout(std_out,msg,'COLL')
334  call wrtout(ab_out,msg,'COLL')
335  write(msg,'(a,i12)')' number of q-points in IBZ                ',Qmesh%nibz
336  call wrtout(std_out,msg,'COLL')
337  call wrtout(ab_out,msg,'COLL')
338  write(msg,'(a,i12)')' number of symmetry operations            ',Cryst%nsym
339  call wrtout(std_out,msg,'COLL')
340  call wrtout(ab_out,msg,'COLL')
341  write(msg,'(a,i12)')' number of k-points in BZ                 ',Kmesh%nbz
342  call wrtout(std_out,msg,'COLL')
343  call wrtout(ab_out,msg,'COLL')
344  write(msg,'(a,i12)')' number of q-points in BZ                 ',Qmesh%nbz
345  call wrtout(std_out,msg,'COLL')
346  call wrtout(ab_out,msg,'COLL')
347  write(msg,'(a,i12)')' number of frequencies for dSigma/dE      ',Sigp%nomegasrd
348  call wrtout(std_out,msg,'COLL')
349  call wrtout(ab_out,msg,'COLL')
350  write(msg,'(a,f12.2)')' frequency step for dSigma/dE [eV]        ',Sigp%deltae*Ha_eV
351  call wrtout(std_out,msg,'COLL')
352  call wrtout(ab_out,msg,'COLL')
353  write(msg,'(a,i12)')' number of omega for Sigma on real axis   ',Sigp%nomegasr
354  call wrtout(std_out,msg,'COLL')
355  call wrtout(ab_out,msg,'COLL')
356  write(msg,'(a,f12.2)')' max omega for Sigma on real axis  [eV]   ',Sigp%maxomega_r*Ha_eV
357  call wrtout(std_out,msg,'COLL')
358  call wrtout(ab_out,msg,'COLL')
359  write(msg,'(a,f12.2)')' zcut for avoiding poles [eV]             ',Sigp%zcut*Ha_eV
360  call wrtout(std_out,msg,'COLL')
361  call wrtout(ab_out,msg,'COLL')
362 
363  if (Sigp%mbpt_sciss>0.1d-4) then
364    write(msg,'(a,f12.2)')' scissor energy [eV]                      ',Sigp%mbpt_sciss*Ha_eV
365    call wrtout(std_out,msg,'COLL')
366    call wrtout(ab_out,msg,'COLL')
367  end if
368 
369  if (mod10==1) then
370    write(msg,'(a,i12)')' number of imaginary frequencies for Sigma',Sigp%nomegasi
371    call wrtout(std_out,msg,'COLL')
372    call wrtout(ab_out,msg,'COLL')
373    write(msg,'(a,f12.2)')' max omega for Sigma on imag axis  [eV]   ',Sigp%omegasimax*Ha_eV
374    call wrtout(std_out,msg,'COLL')
375    call wrtout(ab_out,msg,'COLL')
376  end if
377 
378  if (sigma_needs_w(Sigp)) then
379    write(msg,'(2a)')ch10,' EPSILON^-1 parameters (SCR file):'
380    call wrtout(std_out,msg,'COLL')
381    call wrtout(ab_out,msg,'COLL')
382    !write(std_out,*) titem1(2)(1:79)
383    write(msg,'(a,i12)')' dimension of the eps^-1 matrix on file   ',Er%Hscr%npwe
384    call wrtout(std_out,msg,'COLL')
385    call wrtout(ab_out,msg,'COLL')
386    write(msg,'(a,i12)')' dimension of the eps^-1 matrix used      ',Er%npwe
387    call wrtout(std_out,msg,'COLL')
388    call wrtout(ab_out,msg,'COLL')
389    write(msg,'(a,i12)')' number of plane-waves for wavefunctions  ',Er%Hscr%npwwfn_used
390    call wrtout(std_out,msg,'COLL')
391    call wrtout(ab_out,msg,'COLL')
392    write(msg,'(a,i12)')' number of bands                          ',Er%Hscr%nbnds_used
393    call wrtout(std_out,msg,'COLL')
394    call wrtout(ab_out,msg,'COLL')
395    write(msg,'(a,i12)')' number of q-points in IBZ                ',Qmesh%nibz
396    call wrtout(std_out,msg,'COLL')
397    call wrtout(ab_out,msg,'COLL')
398    write(msg,'(a,i12)')' number of frequencies                    ',Er%nomega
399    call wrtout(std_out,msg,'COLL')
400    call wrtout(ab_out,msg,'COLL')
401    write(msg,'(a,i12)')' number of real frequencies               ',Er%nomega_r
402    call wrtout(std_out,msg,'COLL')
403    call wrtout(ab_out,msg,'COLL')
404    write(msg,'(a,i12)')' number of imag frequencies               ',Er%nomega_i
405    call wrtout(std_out,msg,'COLL')
406    call wrtout(ab_out,msg,'COLL')
407  end if
408 
409  write(msg,'(3a)')ch10,' matrix elements of self-energy operator (all in [eV])',ch10
410  call wrtout(std_out,msg,'COLL')
411  call wrtout(ab_out,msg,'COLL')
412 
413  if (gwcalctyp<10) then
414    write(msg,'(a)')' Perturbative Calculation'
415  else if (gwcalctyp<20) then
416    write(msg,'(a)')' Self-Consistent on Energies only'
417  else
418    write(msg,'(a)')' Self-Consistent on Energies and Wavefunctions'
419  end if
420  call wrtout(std_out,msg,'COLL')
421  call wrtout(ab_out,msg,'COLL')
422 
423 end subroutine write_sigma_header

m_sigma/write_sigma_results [ Functions ]

[ Top ] [ m_sigma ] [ Functions ]

NAME

 write_sigma_results

FUNCTION

  Write the final results of the GW calculation.

INPUTS

  KS_BSt<ebands_t>=Info on the KS band structure energies.
     %eig(mband,nkibz,nsppol)= KS energies
  ikibz= index of the k-point in the array kibz, where GW corrections are calculated
  ikcalc= index of the k-point in the array Sigp%kptgw2bz
  Sigp=sigparams_t datatype
  sr=sigma results datatype

OUTPUT

  (for writing routines, no output) otherwise, should be described

PARENTS

      sigma

CHILDREN

      findqg0,wfd_distribute_bands,wfd_update_bkstab,xmpi_sum

SOURCE

455 subroutine write_sigma_results(ikcalc,ikibz,Sigp,Sr,KS_BSt)
456 
457 
458 !This section has been created automatically by the script Abilint (TD).
459 !Do not modify the following lines by hand.
460 #undef ABI_FUNC
461 #define ABI_FUNC 'write_sigma_results'
462 !End of the abilint section
463 
464  implicit none
465 
466 !Arguments ------------------------------------
467 !scalars
468  integer,intent(in) :: ikcalc,ikibz
469  type(ebands_t),intent(in) :: KS_BSt
470  type(sigparams_t),intent(in) :: Sigp
471  type(sigma_t),intent(in) :: Sr
472 
473 !Local variables-------------------------------
474 !scalars
475  integer :: ib,io,is
476  integer :: gwcalctyp,mod10
477  character(len=500) :: msg
478 !arrays
479  character(len=12) :: tag_spin(2)
480 
481 ! *************************************************************************
482 
483  gwcalctyp=Sigp%gwcalctyp
484  mod10=MOD(Sigp%gwcalctyp,10)
485 
486  !unt_gw  File with GW corrections.
487  !unt_sig Self-energy as a function of frequency.
488  !unt_sgr Derivative wrt omega of the Self-energy.
489  !unt_sgm Sigma on the Matsubara axis.
490 
491  tag_spin=(/'            ','            '/); if (Sr%nsppol==2) tag_spin=(/',  SPIN UP  ',',  SPIN DOWN'/)
492 
493  do is=1,Sr%nsppol
494    write(msg,'(2a,3f8.3,a)')ch10,' k = ',Sigp%kptgw(:,ikcalc),tag_spin(is)
495    call wrtout(std_out,msg,'COLL')
496    call wrtout(ab_out,msg,'COLL')
497 
498    msg = ' Band     E0 <VxcLDA>   SigX SigC(E0)      Z dSigC/dE  Sig(E)    E-E0       E'
499    if (Sr%usepawu/=0) then
500      msg = ' Band     E0 <VxcLDA>   <H_U>  SigX SigC(E0)      Z dSigC/dE  Sig(E)    E-E0       E'
501    end if
502 
503    if (gwcalctyp>=10) then
504      write(msg,'(2a)')&
505 &     ' Band     E_lda   <Vxclda>   E(N-1)  <Hhartree>   SigX  SigC[E(N-1)]',&
506 &     '    Z     dSigC/dE  Sig[E(N)]  DeltaE  E(N)_pert E(N)_diago'
507    end if
508    call wrtout(std_out,msg,'COLL')
509    call wrtout(ab_out,msg,'COLL')
510 
511    write(unt_gw,'(3f10.6)')Sigp%kptgw(:,ikcalc)
512    write(unt_gw,'(i4)')Sigp%maxbnd(ikcalc,is)-Sigp%minbnd(ikcalc,is)+1
513 
514    write(unt_gwdiag,'(3f10.6)')Sigp%kptgw(:,ikcalc)
515    write(unt_gwdiag,'(i4)')Sigp%maxbnd(ikcalc,is)-Sigp%minbnd(ikcalc,is)+1
516 
517    write(unt_sig,'("# k = ",3f10.6)')Sigp%kptgw(:,ikcalc)
518    write(unt_sig,'("# b = ",2i10)')Sigp%minbnd(ikcalc,is),Sigp%maxbnd(ikcalc,is)
519 
520    write(unt_sgr,'("# k = ",3f10.6)')Sigp%kptgw(:,ikcalc)
521    write(unt_sgr,'("# b = ",2i10)')Sigp%minbnd(ikcalc,is),Sigp%maxbnd(ikcalc,is)
522 
523    do ib=Sigp%minbnd(ikcalc,is),Sigp%maxbnd(ikcalc,is)
524      if (gwcalctyp>=10) then
525        call print_Sigma_QPSC(Sr,ikibz,ib,is,KS_BSt,unit=ab_out)
526        call print_Sigma_QPSC(Sr,ikibz,ib,is,KS_BSt,unit=std_out,prtvol=1)
527 
528        write(unt_gwdiag,'(i6,3f9.4)')                                  &
529 &        ib,                                                           &
530 &        Sr%en_qp_diago(ib,ikibz,is)*Ha_eV,                            &
531 &        (Sr%en_qp_diago(ib,ikibz,is) - KS_BSt%eig(ib,ikibz,is))*Ha_eV,&
532 &        zero
533 
534      else
535        ! If not ppmodel, write out also the imaginary part in ab_out
536        SELECT CASE(mod10)
537        CASE(1,2)
538          call print_Sigma_perturbative(Sr,ikibz,ib,is,unit=ab_out,prtvol=1)
539        CASE DEFAULT
540          call print_Sigma_perturbative(Sr,ikibz,ib,is,unit=ab_out)
541        END SELECT
542        call print_Sigma_perturbative(Sr,ikibz,ib,is,unit=std_out,prtvol=1)
543      end if
544 
545      write(unt_gw,'(i6,3f9.4)')          &
546 &      ib,                               &
547 &      REAL (Sr%egw (ib,ikibz,is))*Ha_eV,&
548 &      REAL (Sr%degw(ib,ikibz,is))*Ha_eV,&
549 &      AIMAG(Sr%egw (ib,ikibz,is))*Ha_eV
550    end do !ib
551 
552    if (Sr%e0gap(ikibz,is)**2+Sr%egwgap(ikibz,is)**2+Sr%degwgap(ikibz,is)**2 > tol10) then
553      ! Output the direct gap for each spin
554      ! If all the gaps are zero, this means that it could not be computed in the calling routine
555      write(msg,'(2a,f8.3)')ch10,' E^0_gap       ',Sr%e0gap(ikibz,is)*Ha_eV
556      call wrtout(std_out,msg,'COLL')
557      call wrtout(ab_out,msg,'COLL')
558      write(msg,'(a,f8.3)')      ' E^GW_gap      ',Sr%egwgap(ikibz,is)*Ha_eV
559      call wrtout(std_out,msg,'COLL')
560      call wrtout(ab_out,msg,'COLL')
561      write(msg,'(a,f8.3,a)')    ' DeltaE^GW_gap ',Sr%degwgap(ikibz,is)*Ha_eV,ch10
562      call wrtout(std_out,msg,'COLL')
563      call wrtout(ab_out,msg,'COLL')
564    end if
565    !
566    ! === Output of the spectral function ===
567    do io=1,Sr%nomega_r
568      write(unt_sig,'(100(e12.5,2x))')&
569 &     REAL(Sr%omega_r(io))*Ha_eV,&
570 &     (REAL(Sr%sigxcme(ib,ikibz,io,is))*Ha_eV,&
571 &     AIMAG(Sr%sigxcme(ib,ikibz,io,is))*Ha_eV,&
572 &     gw_spectral_function(Sr,io,ib,ikibz,is),&
573 &     ib=Sigp%minbnd(ikcalc,is),Sigp%maxbnd(ikcalc,is))
574    end do
575    !
576    do ib=Sigp%minbnd(ikcalc,is),Sigp%maxbnd(ikcalc,is)
577      write(unt_sgr,'("# ik, ib",2i5)')ikibz,ib
578      do io=1,Sr%nomega4sd
579        write(unt_sgr,'(100(e12.5,2x))')              &
580 &        REAL (Sr%omega4sd  (ib,ikibz,io,is)) *Ha_eV,&
581 &        REAL (Sr%sigxcme4sd(ib,ikibz,io,is)) *Ha_eV,&
582 &        AIMAG(Sr%sigxcme4sd(ib,ikibz,io,is)) *Ha_eV
583      end do
584    end do
585    !
586    if (mod10==1) then ! For AC, write sigma matrix elements along the imaginary axis
587      do ib=Sigp%minbnd(ikcalc,is),Sigp%maxbnd(ikcalc,is)
588        write(unt_sgm,'("# ik, ib",2i5)')ikibz,ib
589        do io=1,Sr%nomega_i
590          write(unt_sgm,'(3(e12.5,2x))')              &
591 &          AIMAG(Sr%omega_i(io))              *Ha_eV,&
592 &          REAL (Sr%sigxcmesi(ib,ikibz,io,is))*Ha_eV,&
593 &          AIMAG(Sr%sigxcmesi(ib,ikibz,io,is))*Ha_eV
594        end do
595      end do
596    end if
597 
598  end do !is
599 
600 end subroutine write_sigma_results