From bd6fdf87f1087e26d5270d8ceb410f8f91ac2731 Mon Sep 17 00:00:00 2001 From: Pierre-Francois Loos Date: Sun, 23 Jul 2023 23:10:57 +0200 Subject: [PATCH] remove ERI_AO in GW GT --- input/methods | 2 +- src/GT/G0T0eh.f90 | 59 ++++++++++++++++++++++++-------- src/GT/G0T0pp.f90 | 3 +- src/GT/GT.f90 | 13 ++++--- src/GT/GTeh_self_energy_diag.f90 | 2 -- src/GT/UG0T0pp.f90 | 3 +- src/GT/evGTeh.f90 | 3 +- src/GT/evGTpp.f90 | 3 +- src/GT/evUGTpp.f90 | 3 +- src/GW/G0W0.f90 | 21 ++++++------ src/GW/GW.f90 | 10 +++--- src/GW/UG0W0.f90 | 5 ++- src/GW/evGW.f90 | 17 +++++---- src/GW/evUGW.f90 | 3 +- 14 files changed, 83 insertions(+), 64 deletions(-) diff --git a/input/methods b/input/methods index 2305a6a..4ef876d 100644 --- a/input/methods +++ b/input/methods @@ -15,5 +15,5 @@ # G0W0* evGW* qsGW* SRG-qsGW ufG0W0 ufGW F F F F F F # G0T0pp evGTpp qsGTpp G0T0eh evGTeh qsGTeh - F F F F T F + F F F T F F # * unrestricted version available diff --git a/src/GT/G0T0eh.f90 b/src/GT/G0T0eh.f90 index fc7f1a2..65a6c8c 100644 --- a/src/GT/G0T0eh.f90 +++ b/src/GT/G0T0eh.f90 @@ -1,6 +1,6 @@ subroutine G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE,dTDA,doppBSE, & singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF, & - ERI_AO,ERI_MO,dipole_int,PHF,cHF,eHF) + ERI,dipole_int,PHF,cHF,eHF) ! Perform ehG0T0 calculation @@ -34,8 +34,7 @@ subroutine G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE, integer,intent(in) :: nS double precision,intent(in) :: ENuc double precision,intent(in) :: ERHF - double precision,intent(in) :: ERI_AO(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_MO(nBas,nBas,nBas,nBas) + double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int(nBas,nBas,ncart) double precision,intent(in) :: eHF(nBas) double precision,intent(in) :: cHF(nBas,nBas) @@ -43,9 +42,11 @@ subroutine G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE, ! Local variables - logical :: print_W = .true. + logical :: print_T = .true. logical :: dRPA = .false. + logical :: dRPA_W = .true. integer :: ispin + integer :: isp_W double precision :: EcRPA double precision :: EcBSE(nspin) double precision :: EcAC(nspin) @@ -63,6 +64,13 @@ subroutine G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE, double precision,allocatable :: eGT(:) double precision,allocatable :: eGTlin(:) + double precision,allocatable :: KA_sta(:,:) + double precision,allocatable :: KB_sta(:,:) + double precision,allocatable :: OmRPA(:) + double precision,allocatable :: XpY_RPA(:,:) + double precision,allocatable :: XmY_RPA(:,:) + double precision,allocatable :: rho_RPA(:,:,:) + ! Output variables ! Hello world @@ -91,31 +99,54 @@ subroutine G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE, write(*,*) end if -! Spin manifold (triplet for GTeh) - - ispin = 2 - ! Memory allocation allocate(Aph(nS,nS),Bph(nS,nS),Sig(nBas),Z(nBas),Om(nS),XpY(nS,nS),XmY(nS,nS), & rhoL(nBas,nBas,nS),rhoR(nBas,nBas,nS),eGT(nBas),eGTlin(nBas)) +!--------------------------------- +! Compute (singlet) RPA screening +!--------------------------------- + +! allocate(OmRPA(nS),XpY_RPA(nS,nS),XmY_RPA(nS,nS),rho_RPA(nBas,nBas,nS),KA_sta(nS,nS),KB_sta(nS,nS)) + +! isp_W = 1 +! EcRPA = 0d0 + +! call phLR_A(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS,1d0,eHF,ERI,Aph) +! if(.not.TDA_T) call phLR_B(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph) + +! call phLR(TDA_T,nS,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA) +! call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA) + +! call GW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KA_sta) +! call GW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KB_sta) + +! deallocate(OmRPA,XpY_RPA,XmY_RPA,rho_RPA) + !-------------------! ! Compute screening ! !-------------------! - call phLR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,eHF,ERI_MO,Aph) - if(.not.TDA_T) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI_MO,Bph) +! Spin manifold (triplet for GTeh) + + ispin = 2 + + call phLR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,eHF,ERI,Aph) + if(.not.TDA_T) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph) + +! Aph(:,:) = Aph(:,:) + KA_sta(:,:) +! if(.not.TDA_T) Bph(:,:) = Bph(:,:) + KB_sta(:,:) call phLR(TDA_T,nS,Aph,Bph,EcRPA,Om,XpY,XmY) - if(print_W) call print_excitation('RPA@HF ',ispin,nS,Om) + if(print_T) call print_excitation('RPA@HF ',ispin,nS,Om) !--------------------------! ! Compute spectral weights ! !--------------------------! - call GTeh_excitation_density(nBas,nC,nO,nR,nS,ERI_MO,XpY,XmY,rhoL,rhoR) + call GTeh_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,XmY,rhoL,rhoR) !------------------------! ! Compute GW self-energy ! @@ -156,8 +187,8 @@ subroutine G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE, ! Compute the RPA correlation energy based on the G0T0eh quasiparticle energies - call phLR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,eGT,ERI_MO,Aph) - if(.not.TDA_T) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI_MO,Bph) + call phLR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,eGT,ERI,Aph) + if(.not.TDA_T) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph) call phLR(TDA_T,nS,Aph,Bph,EcRPA,Om,XpY,XmY) diff --git a/src/GT/G0T0pp.f90 b/src/GT/G0T0pp.f90 index e1c5c6c..dfdce23 100644 --- a/src/GT/G0T0pp.f90 +++ b/src/GT/G0T0pp.f90 @@ -1,5 +1,5 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,doppBSE,singlet,triplet, & - linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int,PHF,cHF,eHF) + linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_MO,dipole_int,PHF,cHF,eHF) ! Perform one-shot calculation with a T-matrix self-energy (G0T0) @@ -34,7 +34,6 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,dopp double precision,intent(in) :: eHF(nBas) double precision,intent(in) :: cHF(nBas,nBas) double precision,intent(in) :: PHF(nBas,nBas) - double precision,intent(in) :: ERI_AO(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_MO(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int(nBas,nBas,ncart) diff --git a/src/GT/GT.f90 b/src/GT/GT.f90 index 9e2d09f..db7fdca 100644 --- a/src/GT/GT.f90 +++ b/src/GT/GT.f90 @@ -85,11 +85,11 @@ subroutine GT(doG0T0pp,doevGTpp,doqsGTpp,doG0T0eh,doevGTeh,doqsGTeh,unrestricted call wall_time(start_GT) if(unrestricted) then call UG0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,spin_conserved,spin_flip, & - linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb, & + linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_aaaa,ERI_aabb,ERI_bbbb, & dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF) else call G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,doppBSE,singlet,triplet, & - linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI,dipole_int,PHF,cHF,epsHF) + linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI,dipole_int,PHF,cHF,epsHF) end if call wall_time(end_GT) @@ -108,11 +108,11 @@ subroutine GT(doG0T0pp,doevGTpp,doqsGTpp,doG0T0eh,doevGTeh,doqsGTeh,unrestricted call wall_time(start_GT) if(unrestricted) then call evUGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,spin_conserved,spin_flip, & - eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb, & + eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb, & PHF,cHF,epsHF) else call evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,singlet,triplet, & - eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI,dipole_int,PHF,cHF,epsHF) + eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI,dipole_int,PHF,cHF,epsHF) end if call wall_time(end_GT) @@ -157,7 +157,7 @@ subroutine GT(doG0T0pp,doevGTpp,doqsGTpp,doG0T0eh,doevGTeh,doqsGTeh,unrestricted print*,'Unrestricted version of G0T0eh not yet implemented! Sorry.' else call G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE,dTDA,doppBSE,singlet,triplet, & - linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI,dipole_int,PHF,cHF,epsHF) + linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI,dipole_int,PHF,cHF,epsHF) end if call wall_time(end_GT) @@ -178,8 +178,7 @@ subroutine GT(doG0T0pp,doevGTpp,doqsGTpp,doG0T0eh,doevGTeh,doqsGTeh,unrestricted print*,'Unrestricted version of evGTeh not yet implemented! Sorry.' else call evGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE,dTDA,doppBSE, & - singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI,dipole_int, & - PHF,cHF,epsHF) + singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI,dipole_int,PHF,cHF,epsHF) end if call wall_time(end_GT) diff --git a/src/GT/GTeh_self_energy_diag.f90 b/src/GT/GTeh_self_energy_diag.f90 index 43a0350..1fcb612 100644 --- a/src/GT/GTeh_self_energy_diag.f90 +++ b/src/GT/GTeh_self_energy_diag.f90 @@ -47,7 +47,6 @@ subroutine GTeh_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,e,Om,rhoL,rhoR,EcGM,Sig eps = e(p) - e(i) + Om(m) num = rhoL(i,p,m)*rhoR(i,p,m) -! num = 0.5d0*rhoL(p,i,m)*rhoL(p,i,m) Sig(p) = Sig(p) + num*eps/(eps**2 + eta**2) Z(p) = Z(p) - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2 @@ -63,7 +62,6 @@ subroutine GTeh_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,e,Om,rhoL,rhoR,EcGM,Sig eps = e(p) - e(a) - Om(m) num = rhoL(p,a,m)*rhoR(p,a,m) -! num = 0.5d0*rhoL(p,a,m)*rhoL(p,a,m) Sig(p) = Sig(p) + num*eps/(eps**2 + eta**2) Z(p) = Z(p) - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2 diff --git a/src/GT/UG0T0pp.f90 b/src/GT/UG0T0pp.f90 index f6fbb86..a55b099 100644 --- a/src/GT/UG0T0pp.f90 +++ b/src/GT/UG0T0pp.f90 @@ -1,6 +1,6 @@ subroutine UG0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA, & spin_conserved,spin_flip,linearize,eta,regularize,nBas,nC,nO,nV, & - nR,nS,ENuc,EUHF,ERI,ERI_aaaa,ERI_aabb,ERI_bbbb, & + nR,nS,ENuc,EUHF,ERI_aaaa,ERI_aabb,ERI_bbbb, & dipole_int_aa,dipole_int_bb,PHF,cHF,eHF) ! Perform one-shot calculation with a T-matrix self-energy (G0T0) @@ -35,7 +35,6 @@ subroutine UG0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA, & double precision,intent(in) :: eHF(nBas,nspin) double precision,intent(in) :: cHF(nBas,nBas,nspin) double precision,intent(in) :: PHF(nBas,nBas,nspin) - double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) diff --git a/src/GT/evGTeh.f90 b/src/GT/evGTeh.f90 index 9ad7198..554a9c3 100644 --- a/src/GT/evGTeh.f90 +++ b/src/GT/evGTeh.f90 @@ -1,5 +1,5 @@ subroutine evGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE,dTDA,doppBSE, & - singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int,PHF, & + singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_MO,dipole_int,PHF, & cHF,eHF) ! Perform self-consistent eigenvalue-only ehGT calculation @@ -39,7 +39,6 @@ subroutine evGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,d double precision,intent(in) :: PHF(nBas,nBas) double precision,intent(in) :: eHF(nBas) double precision,intent(in) :: cHF(nBas,nBas) - double precision,intent(in) :: ERI_AO(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_MO(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int(nBas,nBas,ncart) diff --git a/src/GT/evGTpp.f90 b/src/GT/evGTpp.f90 index 801ed93..59ed8d8 100644 --- a/src/GT/evGTpp.f90 +++ b/src/GT/evGTpp.f90 @@ -1,5 +1,5 @@ subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA, & - singlet,triplet,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int,PHF,cHF,eHF) + singlet,triplet,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_MO,dipole_int,PHF,cHF,eHF) ! Perform eigenvalue self-consistent calculation with a T-matrix self-energy (evGT) @@ -35,7 +35,6 @@ subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T double precision,intent(in) :: PHF(nBas,nBas) double precision,intent(in) :: eHF(nBas) double precision,intent(in) :: cHF(nBas,nBas) - double precision,intent(in) :: ERI_AO(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_MO(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int(nBas,nBas,ncart) diff --git a/src/GT/evUGTpp.f90 b/src/GT/evUGTpp.f90 index 6b20a5d..2bb9c52 100644 --- a/src/GT/evUGTpp.f90 +++ b/src/GT/evUGTpp.f90 @@ -1,6 +1,6 @@ subroutine evUGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE, & TDA_T,TDA,dBSE,dTDA,spin_conserved,spin_flip,& - eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,ERI_AO,ERI_aaaa, & + eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,ERI_aaaa, & ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF) ! Perform one-shot calculation with a T-matrix self-energy (G0T0) @@ -36,7 +36,6 @@ subroutine evUGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE, & double precision,intent(in) :: eHF(nBas,nspin) double precision,intent(in) :: cHF(nBas,nBas,nspin) double precision,intent(in) :: PHF(nBas,nBas,nspin) - double precision,intent(in) :: ERI_AO(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) diff --git a/src/GW/G0W0.f90 b/src/GW/G0W0.f90 index 5ae8a1b..29a026f 100644 --- a/src/GW/G0W0.f90 +++ b/src/GW/G0W0.f90 @@ -1,5 +1,5 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE,singlet,triplet, & - linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int,PHF,cHF,eHF) + linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,dipole_int,PHF,cHF,eHF) ! Perform G0W0 calculation @@ -33,8 +33,7 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dT integer,intent(in) :: nS double precision,intent(in) :: ENuc double precision,intent(in) :: ERHF - double precision,intent(in) :: ERI_AO(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_MO(nBas,nBas,nBas,nBas) + double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int(nBas,nBas,ncart) double precision,intent(in) :: eHF(nBas) double precision,intent(in) :: cHF(nBas,nBas) @@ -102,8 +101,8 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dT ! Compute screening ! !-------------------! - call phLR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,eHF,ERI_MO,Aph) - if(.not.TDA_W) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI_MO,Bph) + call phLR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,eHF,ERI,Aph) + if(.not.TDA_W) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph) call phLR(TDA_W,nS,Aph,Bph,EcRPA,Om,XpY,XmY) @@ -113,7 +112,7 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dT ! Compute spectral weights ! !--------------------------! - call GW_excitation_density(nBas,nC,nO,nR,nS,ERI_MO,XpY,rho) + call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho) !------------------------! ! Compute GW self-energy ! @@ -156,8 +155,8 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dT ! Compute the RPA correlation energy - call phLR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,eGW,ERI_MO,Aph) - if(.not.TDA_W) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI_MO,Bph) + call phLR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,eGW,ERI,Aph) + if(.not.TDA_W) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph) call phLR(TDA_W,nS,Aph,Bph,EcRPA,Om,XpY,XmY) @@ -175,7 +174,7 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dT if(dophBSE) then - call GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int,eHF,eGW,EcBSE) + call GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGW,EcBSE) if(exchange_kernel) then @@ -210,7 +209,7 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dT end if - call GW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,dophBSE,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,eHF,eGW,EcBSE) + call GW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,dophBSE,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,eHF,eGW,EcBSE) write(*,*) write(*,*)'-------------------------------------------------------------------------------' @@ -227,7 +226,7 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dT if(doppBSE) then - call GW_ppBSE(TDA_W,TDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int,eHF,eGW,EcBSE) + call GW_ppBSE(TDA_W,TDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGW,EcBSE) write(*,*) write(*,*)'-------------------------------------------------------------------------------' diff --git a/src/GW/GW.f90 b/src/GW/GW.f90 index 8e92d5e..ff5692c 100644 --- a/src/GW/GW.f90 +++ b/src/GW/GW.f90 @@ -83,12 +83,12 @@ subroutine GW(doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,unrestricted,maxSCF call wall_time(start_GW) if(unrestricted) then - call UG0W0(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip, & - linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,S,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb, & + call UG0W0(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip, & + linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,S,ERI_aaaa,ERI_aabb,ERI_bbbb, & dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF) else call G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE,singlet,triplet, & - linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI,dipole_int,PHF,cHF,epsHF) + linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI,dipole_int,PHF,cHF,epsHF) end if call wall_time(end_GW) @@ -107,11 +107,11 @@ subroutine GW(doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,unrestricted,maxSCF call wall_time(start_GW) if(unrestricted) then call evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip, & - eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,S,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb, & + eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,S,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb, & PHF,cHF,epsHF) else call evGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, & - singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI,dipole_int,PHF,cHF,epsHF) + singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI,dipole_int,PHF,cHF,epsHF) end if call wall_time(end_GW) diff --git a/src/GW/UG0W0.f90 b/src/GW/UG0W0.f90 index 79600f1..9ae32e4 100644 --- a/src/GW/UG0W0.f90 +++ b/src/GW/UG0W0.f90 @@ -1,5 +1,5 @@ -subroutine UG0W0(doACFDT,exchange_kernel,doXBS,BSE,TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip, & - linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,S,ERI,ERI_aaaa,ERI_aabb,ERI_bbbb, & +subroutine UG0W0(doACFDT,exchange_kernel,doXBS,BSE,TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip, & + linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,S,ERI_aaaa,ERI_aabb,ERI_bbbb, & dipole_int_aa,dipole_int_bb,PHF,cHF,eHF) ! Perform unrestricted G0W0 calculation @@ -33,7 +33,6 @@ subroutine UG0W0(doACFDT,exchange_kernel,doXBS,BSE,TDA_W,TDA,dBSE,dTDA,spin_cons double precision,intent(in) :: ENuc double precision,intent(in) :: EUHF double precision,intent(in) :: S(nBas,nBas) - double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) diff --git a/src/GW/evGW.f90 b/src/GW/evGW.f90 index 9f17631..a763a90 100644 --- a/src/GW/evGW.f90 +++ b/src/GW/evGW.f90 @@ -1,5 +1,5 @@ subroutine evGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, & - singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int,PHF, & + singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,dipole_int,PHF, & cHF,eHF) ! Perform self-consistent eigenvalue-only GW calculation @@ -39,8 +39,7 @@ subroutine evGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dop double precision,intent(in) :: PHF(nBas,nBas) double precision,intent(in) :: eHF(nBas) double precision,intent(in) :: cHF(nBas,nBas) - double precision,intent(in) :: ERI_AO(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_MO(nBas,nBas,nBas,nBas) + double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int(nBas,nBas,ncart) ! Local variables @@ -122,14 +121,14 @@ subroutine evGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dop ! Compute screening - call phLR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,eGW,ERI_MO,Aph) - if(.not.TDA_W) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI_MO,Bph) + call phLR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,eGW,ERI,Aph) + if(.not.TDA_W) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph) call phLR(TDA_W,nS,Aph,Bph,EcRPA,Om,XpY,XmY) ! Compute spectral weights - call GW_excitation_density(nBas,nC,nO,nR,nS,ERI_MO,XpY,rho) + call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho) ! Compute correlation part of the self-energy @@ -226,7 +225,7 @@ subroutine evGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dop if(dophBSE) then - call GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int,eGW,eGW,EcBSE) + call GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eGW,eGW,EcBSE) if(exchange_kernel) then @@ -260,7 +259,7 @@ subroutine evGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dop end if - call GW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,dophBSE,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,eGW,eGW,EcBSE) + call GW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,dophBSE,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,eGW,eGW,EcBSE) write(*,*) write(*,*)'-------------------------------------------------------------------------------' @@ -277,7 +276,7 @@ subroutine evGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dop if(doppBSE) then - call GW_ppBSE(TDA_W,TDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int,eHF,eGW,EcBSE) + call GW_ppBSE(TDA_W,TDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGW,EcBSE) write(*,*) write(*,*)'-------------------------------------------------------------------------------' diff --git a/src/GW/evUGW.f90 b/src/GW/evUGW.f90 index e36bf0f..0030bfc 100644 --- a/src/GW/evUGW.f90 +++ b/src/GW/evUGW.f90 @@ -1,6 +1,6 @@ subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_W,TDA, & dBSE,dTDA,spin_conserved,spin_flip,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc, & - EUHF,S,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF) + EUHF,S,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF) ! Perform self-consistent eigenvalue-only GW calculation @@ -38,7 +38,6 @@ subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_W, double precision,intent(in) :: eHF(nBas,nspin) double precision,intent(in) :: cHF(nBas,nBas,nspin) double precision,intent(in) :: S(nBas,nBas) - double precision,intent(in) :: ERI_AO(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas)