TABLE OF CONTENTS
ABINIT/spectral_function [ Functions ]
NAME
spectral_function
FUNCTION
Print the spectral function computed from Green function in real frequency
COPYRIGHT
Copyright (C) 1999-2018 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
cryst_struc <type(crystal_t)>=crystal structure data green <type(green_type)>= green function data hu <type(hu_type)>= datatype of type hu paw_dmft = data for self-consistent LDA+DMFT calculations. pawang <type(pawang)>=paw angular mesh and related data self <type(self_type)>= variables related to self-energy prtopt= option for printing
OUTPUT
paw_dmft = data for self-consistent LDA+DMFT calculations.
NOTES
PARENTS
dmft_solve
CHILDREN
compute_green,copy_green,copy_matlu,dc_self,destroy_green,destroy_self dyson,hubbard_one,init_green,initialize_self,ldau_self,print_green rw_self,wrtout
SOURCE
40 #if defined HAVE_CONFIG_H 41 #include "config.h" 42 #endif 43 44 #include "abi_common.h" 45 46 subroutine spectral_function(cryst_struc,green,hu,paw_dmft,& 47 & pawang,pawtab,self_old,prtopt) 48 49 use defs_basis 50 use defs_abitypes 51 use m_errors 52 use m_profiling_abi 53 54 use m_crystal, only : crystal_t 55 use m_green, only : init_green,green_type,print_green,copy_green,compute_green,destroy_green 56 use m_matlu, only : copy_matlu 57 use m_paw_dmft, only : paw_dmft_type 58 use m_hu, only : hu_type 59 use m_self, only : self_type,initialize_self,dc_self,destroy_self,rw_self 60 use m_energy, only : energy_type 61 use m_pawang, only : pawang_type 62 use m_pawtab, only : pawtab_type 63 64 !This section has been created automatically by the script Abilint (TD). 65 !Do not modify the following lines by hand. 66 #undef ABI_FUNC 67 #define ABI_FUNC 'spectral_function' 68 use interfaces_14_hidewrite 69 use interfaces_68_dmft, except_this_one => spectral_function 70 !End of the abilint section 71 72 implicit none 73 74 !Arguments ------------------------------------ 75 !scalars 76 ! type(pawang_type), intent(in) :: pawang 77 type(crystal_t),intent(in) :: cryst_struc 78 type(green_type), intent(in) :: green 79 type(hu_type),intent(inout) :: hu(cryst_struc%ntypat) 80 !type(MPI_type), intent(inout) :: mpi_enreg 81 type(pawang_type), intent(in) :: pawang 82 type(pawtab_type),intent(in) :: pawtab(cryst_struc%ntypat) 83 type(self_type), intent(inout) :: self_old 84 type(paw_dmft_type), intent(inout) :: paw_dmft 85 integer, intent(in) :: prtopt 86 87 !Local variables ------------------------------ 88 character(len=500) :: message 89 type(green_type) :: greenr 90 type(green_type) :: weissr 91 type(self_type) :: selfr 92 !scalars 93 !************************************************************************ 94 !character(len=500) :: message 95 96 ! opt_oper_ksloc=3 to be able to compute spectral function. 97 call init_green(greenr,paw_dmft,opt_oper_ksloc=3,wtype="real") 98 call init_green(weissr,paw_dmft,wtype="real") 99 call copy_matlu(green%occup%matlu,greenr%occup%matlu,paw_dmft%natom) 100 call initialize_self(selfr,paw_dmft,wtype="real") 101 !======================================================================= 102 !== Solve impurity model with green function for real frequency 103 !======================================================================= 104 write(message,'(2a,i3,13x,a)') ch10,' === Write Spectral function' 105 call wrtout(std_out,message,'COLL') 106 if(abs(paw_dmft%dmft_solv)==1) then 107 108 ! == LDA+U for test 109 ! ------------------- 110 call ldau_self(cryst_struc,greenr,paw_dmft,& 111 & pawtab,selfr,opt_ldau=1,prtopt=prtopt) 112 else if(abs(paw_dmft%dmft_solv)==2) then 113 114 ! == Hubbard One 115 ! ------------------- 116 call hubbard_one(cryst_struc,greenr,hu,paw_dmft,& 117 & pawang,prtopt,self_old%hdc,weissr) 118 119 else if(abs(paw_dmft%dmft_solv)==4) then 120 121 ! == Nothing 122 ! ------------------- 123 message = "spectral_function: This section of code is disabled!" 124 MSG_ERROR(message) 125 call copy_green(weissr,greenr,opt_tw=1) 126 127 else if(abs(paw_dmft%dmft_solv)>=5) then 128 129 ! == Nothing 130 ! ------------------- 131 MSG_ERROR("Stopping before copy_green") 132 call copy_green(weissr,greenr,opt_tw=1) 133 134 else if(abs(paw_dmft%dmft_solv)==0) then 135 136 ! == Nothing 137 ! ------------------- 138 ! weiss%occup%has_operks=0 -> only local part is duplicated 139 call copy_green(weissr,greenr,opt_tw=2) 140 end if 141 142 !======================================================================= 143 !== Integrate green function and printout occupations 144 !For dmft_solv=-1,0,or 1 , the green function was not computed: it 145 !cannot be integrated 146 !======================================================================= 147 call dc_self(green%charge_matlu_solver,cryst_struc,hu,selfr,paw_dmft%dmft_dc,prtopt) 148 if(abs(paw_dmft%dmft_solv)/=1.and.paw_dmft%dmft_solv/=0) then 149 call dyson(greenr,paw_dmft,selfr,weissr,opt_weissself=2) 150 end if 151 call compute_green(cryst_struc,greenr,paw_dmft,pawang,1,selfr,opt_self=1) 152 call print_green("realw",greenr,4,paw_dmft,pawprtvol=3) 153 call rw_self(selfr,paw_dmft,prtopt=2,opt_rw=2) 154 155 if(abs(prtopt)>0) then 156 end if 157 call destroy_self(selfr) 158 call destroy_green(weissr) 159 call destroy_green(greenr) 160 161 end subroutine spectral_function