From c2fc694fdb48136da448cc49f0e312cce11f6419 Mon Sep 17 00:00:00 2001 From: pfloos Date: Mon, 29 Jan 2024 15:23:34 +0100 Subject: [PATCH] GW spectal function --- src/GW/GW_plot_self_energy.f90 | 20 +++++++++++--------- 1 file changed, 11 insertions(+), 9 deletions(-) diff --git a/src/GW/GW_plot_self_energy.f90 b/src/GW/GW_plot_self_energy.f90 index 5c60440..0978c39 100644 --- a/src/GW/GW_plot_self_energy.f90 +++ b/src/GW/GW_plot_self_energy.f90 @@ -24,9 +24,9 @@ subroutine GW_plot_self_energy(nBas,nC,nO,nV,nR,nS,eHF,eGW,Om,rho) integer :: p,g integer :: nGrid double precision :: wmin,wmax,dw - double precision,external :: GW_SigC,GW_dSigC + double precision,external :: GW_SigC,GW_ImSigC,GW_dSigC double precision,allocatable :: w(:) - double precision,allocatable :: SigC(:,:) + double precision,allocatable :: ReSigC(:,:),ImSigC(:,:) double precision,allocatable :: Z(:,:) double precision,allocatable :: S(:,:) @@ -37,12 +37,13 @@ subroutine GW_plot_self_energy(nBas,nC,nO,nV,nR,nS,eHF,eGW,Om,rho) ! Construct grid nGrid = 5000 - allocate(w(nGrid),SigC(nBas,nGrid),Z(nBas,nGrid),S(nBas,nGrid)) + allocate(w(nGrid),ReSigC(nBas,nGrid),ImSigC(nBas,nGrid),Z(nBas,nGrid),S(nBas,nGrid)) ! Initialize - SigC(:,:) = 0d0 - Z(:,:) = 0d0 + ReSigC(:,:) = 0d0 + ImSigC(:,:) = 0d0 + Z(:,:) = 0d0 ! Minimum and maximum frequency values @@ -59,8 +60,9 @@ subroutine GW_plot_self_energy(nBas,nC,nO,nV,nR,nS,eHF,eGW,Om,rho) do g=1,nGrid do p=nC+1,nBas-nR - SigC(p,g) = GW_SigC(p,w(g),eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho) - Z(p,g) = GW_dSigC(p,w(g),eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho) + ReSigC(p,g) = GW_SigC(p,w(g),eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho) + ImSigC(p,g) = GW_ImSigC(p,w(g),eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho) + Z(p,g) = GW_dSigC(p,w(g),eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho) end do end do @@ -71,7 +73,7 @@ subroutine GW_plot_self_energy(nBas,nC,nO,nV,nR,nS,eHF,eGW,Om,rho) do g=1,nGrid do p=nC+1,nBas-nR - S(p,g) = eta/((w(g) - eHF(p) - SigC(p,g))**2 + eta**2) + S(p,g) = abs(ImSigC(p,g))/((w(g) - eHF(p) - ReSigC(p,g))**2 + ImSigC(p,g)**2) end do end do @@ -85,7 +87,7 @@ subroutine GW_plot_self_energy(nBas,nC,nO,nV,nR,nS,eHF,eGW,Om,rho) open(unit=11 ,file='GW_A.dat') do g=1,nGrid - write(8 ,*) w(g)*HaToeV,(SigC(p,g)*HaToeV,p=nC+1,nBas-nR) + write(8 ,*) w(g)*HaToeV,(ReSigC(p,g)*HaToeV,p=nC+1,nBas-nR) write(9 ,*) w(g)*HaToeV,((w(g)-eHF(p))*HaToeV,p=nC+1,nBas-nR) write(10,*) w(g)*HaToeV,(Z(p,g),p=nC+1,nBas-nR) write(11,*) w(g)*HaToeV,(S(p,g),p=nC+1,nBas-nR)