From ba7a2349afd37dbebca4961639f60b6f7300735c Mon Sep 17 00:00:00 2001 From: Pierre-Francois Loos Date: Mon, 15 Feb 2021 22:15:30 +0100 Subject: [PATCH] fix QP graph --- src/MBPT/G0W0.f90 | 2 +- src/MBPT/QP_graph.f90 | 16 +++++++++------- src/MBPT/UG0W0.f90 | 4 ++-- src/MBPT/unrestricted_QP_graph.f90 | 16 +++++++++------- 4 files changed, 21 insertions(+), 17 deletions(-) diff --git a/src/MBPT/G0W0.f90 b/src/MBPT/G0W0.f90 index d2150b1..d2c7fa7 100644 --- a/src/MBPT/G0W0.f90 +++ b/src/MBPT/G0W0.f90 @@ -155,7 +155,7 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,COHSEX,SOSEX,BSE,TDA_W,TDA, & write(*,*) ' *** Quasiparticle energies obtained by root search (experimental) *** ' write(*,*) - call QP_graph(nBas,nC,nO,nV,nR,nS,eta,eHF,OmRPA,rho_RPA,eG0W0lin,eG0W0) + call QP_graph(nBas,nC,nO,nV,nR,nS,eta,eHF,SigX,Vxc,OmRPA,rho_RPA,eG0W0lin,eG0W0) ! Find all the roots of the QP equation if necessary diff --git a/src/MBPT/QP_graph.f90 b/src/MBPT/QP_graph.f90 index bc6ae48..69628e8 100644 --- a/src/MBPT/QP_graph.f90 +++ b/src/MBPT/QP_graph.f90 @@ -1,4 +1,4 @@ -subroutine QP_graph(nBas,nC,nO,nV,nR,nS,eta,eHF,Omega,rho,eGWlin,eGW) +subroutine QP_graph(nBas,nC,nO,nV,nR,nS,eta,eHF,SigX,Vxc,Omega,rho,eGWlin,eGW) ! Compute the graphical solution of the QP equation @@ -15,6 +15,8 @@ subroutine QP_graph(nBas,nC,nO,nV,nR,nS,eta,eHF,Omega,rho,eGWlin,eGW) integer,intent(in) :: nS double precision,intent(in) :: eta double precision,intent(in) :: eHF(nBas) + double precision,intent(in) :: SigX(nBas) + double precision,intent(in) :: Vxc(nBas) double precision,intent(in) :: Omega(nS) double precision,intent(in) :: rho(nBas,nBas,nS) @@ -27,7 +29,7 @@ subroutine QP_graph(nBas,nC,nO,nV,nR,nS,eta,eHF,Omega,rho,eGWlin,eGW) integer,parameter :: maxIt = 64 double precision,parameter :: thresh = 1d-6 double precision,external :: SigmaC,dSigmaC - double precision :: sig,dsig + double precision :: sigC,dsigC double precision :: f,df double precision :: w @@ -52,14 +54,14 @@ subroutine QP_graph(nBas,nC,nO,nV,nR,nS,eta,eHF,Omega,rho,eGWlin,eGW) nIt = nIt + 1 - sig = SigmaC(p,w,eta,nBas,nC,nO,nV,nR,nS,eHF,Omega,rho) - dsig = dSigmaC(p,w,eta,nBas,nC,nO,nV,nR,nS,eHF,Omega,rho) - f = w - eHF(p) - sig - df = 1d0 - dsig + sigC = SigmaC(p,w,eta,nBas,nC,nO,nV,nR,nS,eHF,Omega,rho) + dsigC = dSigmaC(p,w,eta,nBas,nC,nO,nV,nR,nS,eHF,Omega,rho) + f = w - eHF(p) - SigX(p) - sigC + Vxc(p) + df = 1d0 - dsigC w = w - f/df - write(*,'(A3,I3,A1,1X,3F15.9)') 'It.',nIt,':',w*HaToeV,f,sig + write(*,'(A3,I3,A1,1X,3F15.9)') 'It.',nIt,':',w*HaToeV,f,sigC end do diff --git a/src/MBPT/UG0W0.f90 b/src/MBPT/UG0W0.f90 index 5c876cc..5fd0131 100644 --- a/src/MBPT/UG0W0.f90 +++ b/src/MBPT/UG0W0.f90 @@ -165,8 +165,8 @@ subroutine UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,ev ! Find graphical solution of the QP equation do is=1,nspin - call unrestricted_QP_graph(nBas,nC(is),nO(is),nV(is),nR(is),nS_sc,eta,eHF(:,is),OmRPA, & - rho_RPA(:,:,:,is),eGWlin(:,is),eGW(:,is)) + call unrestricted_QP_graph(nBas,nC(is),nO(is),nV(is),nR(is),nS_sc,eta,eHF(:,is),SigX(:,is),Vxc(:,is), & + OmRPA,rho_RPA(:,:,:,is),eGWlin(:,is),eGW(:,is)) end do end if diff --git a/src/MBPT/unrestricted_QP_graph.f90 b/src/MBPT/unrestricted_QP_graph.f90 index 96585a9..14b7653 100644 --- a/src/MBPT/unrestricted_QP_graph.f90 +++ b/src/MBPT/unrestricted_QP_graph.f90 @@ -1,4 +1,4 @@ -subroutine unrestricted_QP_graph(nBas,nC,nO,nV,nR,nS,eta,eHF,Omega,rho,eGWlin,eGW) +subroutine unrestricted_QP_graph(nBas,nC,nO,nV,nR,nS,eta,eHF,SigX,Vxc,Omega,rho,eGWlin,eGW) ! Compute the graphical solution of the QP equation @@ -15,6 +15,8 @@ subroutine unrestricted_QP_graph(nBas,nC,nO,nV,nR,nS,eta,eHF,Omega,rho,eGWlin,eG integer,intent(in) :: nS double precision,intent(in) :: eta double precision,intent(in) :: eHF(nBas) + double precision,intent(in) :: SigX(nBas) + double precision,intent(in) :: Vxc(nBas) double precision,intent(in) :: Omega(nS) double precision,intent(in) :: rho(nBas,nBas,nS,nspin) @@ -27,7 +29,7 @@ subroutine unrestricted_QP_graph(nBas,nC,nO,nV,nR,nS,eta,eHF,Omega,rho,eGWlin,eG integer,parameter :: maxIt = 10 double precision,parameter :: thresh = 1d-6 double precision,external :: USigmaC,dUSigmaC - double precision :: sig,dsig + double precision :: sigC,dsigC double precision :: f,df double precision :: w @@ -52,14 +54,14 @@ subroutine unrestricted_QP_graph(nBas,nC,nO,nV,nR,nS,eta,eHF,Omega,rho,eGWlin,eG nIt = nIt + 1 - sig = USigmaC(p,w,eta,nBas,nC,nO,nV,nR,nS,eHF,Omega,rho) - dsig = dUSigmaC(p,w,eta,nBas,nC,nO,nV,nR,nS,eHF,Omega,rho) - f = w - eHF(p) - sig - df = 1d0 - dsig + sigC = USigmaC(p,w,eta,nBas,nC,nO,nV,nR,nS,eHF,Omega,rho) + dsigC = dUSigmaC(p,w,eta,nBas,nC,nO,nV,nR,nS,eHF,Omega,rho) + f = w - eHF(p) - SigX(p) + Vxc(p) - sigC + df = 1d0 - dsigC w = w - f/df - write(*,'(A3,I3,A1,1X,3F15.9)') 'It.',nIt,':',w*HaToeV,f,sig + write(*,'(A3,I3,A1,1X,3F15.9)') 'It.',nIt,':',w*HaToeV,f,sigC end do