diff --git a/input/dft b/input/dft index 3220f9d..8059451 100644 --- a/input/dft +++ b/input/dft @@ -1,12 +1,12 @@ # Restricted or unrestricted KS calculation - UKS + eDFT-UKS # exchange rung: # Hartree = 0: H # LDA = 1: S51,CC-S51 # GGA = 2: B88,G96,PBE # MGGA = 3: # Hybrid = 4: HF,B3,PBE - 1 S51 + 4 HF # correlation rung: # Hartree = 0: H # LDA = 1: PW92,VWN3,VWN5,eVWN5 @@ -17,18 +17,18 @@ # quadrature grid SG-n 1 # Number of states in ensemble (nEns) - 1 + 2 # occupation numbers - 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 - 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 + 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 + 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 - 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 - 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 + 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 + 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 - 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 - 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 + 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 + 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 # Ensemble weights: wEns(1),...,wEns(nEns-1) - 0.0 0.0 + 0.5 0.0 # N-centered? F # Parameters for CC weight-dependent exchange functional diff --git a/input/methods b/input/methods index 7fa0bc0..098169e 100644 --- a/input/methods +++ b/input/methods @@ -11,7 +11,7 @@ # RPA* RPAx* ppRPA F F F # G0F2* evGF2* qsGF2* G0F3 evGF3 - F T F F F + F F T F F # G0W0* evGW* qsGW* F F F # G0T0 evGT qsGT diff --git a/src/GF/print_qsGF2.f90 b/src/GF/print_qsGF2.f90 index 5e3faf9..33a6e9b 100644 --- a/src/GF/print_qsGF2.f90 +++ b/src/GF/print_qsGF2.f90 @@ -69,7 +69,7 @@ subroutine print_qsGF2(nBas,nO,nSCF,Conv,thresh,eHF,eGF2,c,P,T,V,J,K,F,SigC,Z, & write(*,'(2X,A30,F15.6,A3)') 'qsGF2 LUMO energy:',eGF2(LUMO)*HaToeV,' eV' write(*,'(2X,A30,F15.6,A3)') 'qsGF2 HOMO-LUMO gap :',Gap*HaToeV,' eV' write(*,*)'-------------------------------------------' - write(*,'(2X,A30,F15.6,A3)') ' qsGF2 total energy:',ENuc + EqsGF2 + Ec,' au' + write(*,'(2X,A30,F15.6,A3)') ' qsGF2 total energy:',ENuc + EqsGF2,' au' write(*,'(2X,A30,F15.6,A3)') ' qsGF2 exchange energy:',Ex,' au' write(*,'(2X,A30,F15.6,A3)') ' qsGF2 correlation energy:',Ec,' au' write(*,*)'-------------------------------------------' @@ -87,14 +87,14 @@ subroutine print_qsGF2(nBas,nO,nSCF,Conv,thresh,eHF,eGF2,c,P,T,V,J,K,F,SigC,Z, & write(*,'(A32,1X,F16.10,A3)') ' Kinetic energy: ',ET,' au' write(*,'(A32,1X,F16.10,A3)') ' Potential energy: ',EV,' au' write(*,'(A50)') '---------------------------------------' - write(*,'(A32,1X,F16.10,A3)') ' Two-electron energy: ',EJ + Ex,' au' + write(*,'(A32,1X,F16.10,A3)') ' Two-electron energy: ',EJ + Ex + Ec,' au' write(*,'(A32,1X,F16.10,A3)') ' Hartree energy: ',EJ,' au' write(*,'(A32,1X,F16.10,A3)') ' Exchange energy: ',Ex,' au' write(*,'(A32,1X,F16.10,A3)') ' Correlation energy: ',Ec,' au' write(*,'(A50)') '---------------------------------------' - write(*,'(A32,1X,F16.10,A3)') ' Electronic energy: ',EqsGF2 + Ec,' au' + write(*,'(A32,1X,F16.10,A3)') ' Electronic energy: ',EqsGF2,' au' write(*,'(A32,1X,F16.10,A3)') ' Nuclear repulsion: ',ENuc,' au' - write(*,'(A32,1X,F16.10,A3)') ' qsGF2 energy: ',ENuc + EqsGF2 + Ec,' au' + write(*,'(A32,1X,F16.10,A3)') ' qsGF2 energy: ',ENuc + EqsGF2,' au' write(*,'(A50)') '---------------------------------------' write(*,'(A35)') ' Dipole moment (Debye) ' write(*,'(10X,4A10)') 'X','Y','Z','Tot.' diff --git a/src/GF/print_qsUGF2.f90 b/src/GF/print_qsUGF2.f90 index 4db8f23..8c83812 100644 --- a/src/GF/print_qsUGF2.f90 +++ b/src/GF/print_qsUGF2.f90 @@ -99,7 +99,7 @@ subroutine print_qsUGF2(nBas,nO,nSCF,Conv,thresh,eHF,eGF2,cGF2,PGF2,Ov,T,V,J,K, write(*,'(2X,A30,F15.6,A3)') 'qsUGF2 HOMO-LUMO gap :',(minval(LUMO(:))-maxval(HOMO(:)))*HaToeV,' eV' write(*,*)'-------------------------------------------------------------------------------& -------------------------------------------------' - write(*,'(2X,A30,F15.6,A3)') ' qsUGF2 total energy:',ENuc + EqsGF2 + sum(Ec(:)),' au' + write(*,'(2X,A30,F15.6,A3)') ' qsUGF2 total energy:',ENuc + EqsGF2,' au' write(*,'(2X,A30,F15.6,A3)') ' qsUGF2 exchange energy:',sum(Ex(:)),' au' write(*,'(2X,A30,F15.6,A3)') ' qsUGF2 correlation energy:',sum(Ec(:)),' au' write(*,*)'-------------------------------------------------------------------------------& @@ -126,10 +126,10 @@ subroutine print_qsUGF2(nBas,nO,nSCF,Conv,thresh,eHF,eGF2,cGF2,PGF2,Ov,T,V,J,K, write(*,'(A40,1X,F16.10,A3)') ' Potential a energy: ',EV(1),' au' write(*,'(A40,1X,F16.10,A3)') ' Potential b energy: ',EV(2),' au' write(*,'(A60)') '-------------------------------------------------' - write(*,'(A40,1X,F16.10,A3)') ' Two-electron energy: ',sum(EJ(:)) + sum(Ex(:)),' au' - write(*,'(A40,1X,F16.10,A3)') ' Two-electron aa energy: ',EJ(1) + Ex(1),' au' - write(*,'(A40,1X,F16.10,A3)') ' Two-electron ab energy: ',EJ(2),' au' - write(*,'(A40,1X,F16.10,A3)') ' Two-electron bb energy: ',EJ(3) + Ex(2),' au' + write(*,'(A40,1X,F16.10,A3)') ' Two-electron energy: ',sum(EJ(:)) + sum(Ex(:)) + sum(Ec(:)),' au' + write(*,'(A40,1X,F16.10,A3)') ' Two-electron aa energy: ',EJ(1) + Ex(1) + Ec(1),' au' + write(*,'(A40,1X,F16.10,A3)') ' Two-electron ab energy: ',EJ(2) + Ec(2),' au' + write(*,'(A40,1X,F16.10,A3)') ' Two-electron bb energy: ',EJ(3) + Ex(2) + Ec(3),' au' write(*,*) write(*,'(A40,1X,F16.10,A3)') ' Hartree energy: ',sum(EJ(:)),' au' write(*,'(A40,1X,F16.10,A3)') ' Hartree aa energy: ',EJ(1),' au' diff --git a/src/GF/qsGF2.f90 b/src/GF/qsGF2.f90 index 2f39f70..6e8988b 100644 --- a/src/GF/qsGF2.f90 +++ b/src/GF/qsGF2.f90 @@ -197,16 +197,17 @@ subroutine qsGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,evDyn,singlet,triplet, ! Exchange energy - Ex = -0.25d0*trace_matrix(nBas,matmul(P,K)) - - ! Total energy - - EqsGF2 = ET + EV + EJ + Ex + Ex = 0.25d0*trace_matrix(nBas,matmul(P,K)) ! Correlation energy call MP2(nBas,nC,nO,nV,nR,ERI_MO,ENuc,EqsGF2,eGF2,Ec) + ! Total energy + + EqsGF2 = ET + EV + EJ + Ex + Ec + + !------------------------------------------------------------------------ ! Print results !------------------------------------------------------------------------ @@ -248,8 +249,8 @@ subroutine qsGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,evDyn,singlet,triplet, write(*,*)'-------------------------------------------------------------------------------' write(*,'(2X,A50,F20.10)') 'Tr@BSE@qsGF2 correlation energy (singlet) =',EcBSE(1) write(*,'(2X,A50,F20.10)') 'Tr@BSE@qsGF2 correlation energy (triplet) =',EcBSE(2) - write(*,'(2X,A50,F20.10)') 'Tr@BSE@qsGF2 correlation energy =',EcBSE(1) + EcBSE(2) - write(*,'(2X,A50,F20.10)') 'Tr@BSE@qsGF2 total energy =',ENuc + EqsGF2 + EcBSE(1) + EcBSE(2) + write(*,'(2X,A50,F20.10)') 'Tr@BSE@qsGF2 correlation energy =',sum(EcBSE(:)) + write(*,'(2X,A50,F20.10)') 'Tr@BSE@qsGF2 total energy =',ENuc + EqsGF2 + sum(EcBSE(:)) write(*,*)'-------------------------------------------------------------------------------' write(*,*) diff --git a/src/GF/qsUGF2.f90 b/src/GF/qsUGF2.f90 index 70a8e68..5c93e01 100644 --- a/src/GF/qsUGF2.f90 +++ b/src/GF/qsUGF2.f90 @@ -64,7 +64,7 @@ subroutine qsUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,evDyn,spin_conserved, double precision :: EV(nspin) double precision :: EJ(nsp) double precision :: Ex(nspin) - double precision :: Ec(nspin) + double precision :: Ec(nsp) double precision :: EqsGF2 double precision :: EcBSE(nspin) double precision :: EcAC(nspin) @@ -268,14 +268,14 @@ subroutine qsUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,evDyn,spin_conserved, Ex(is) = 0.5d0*trace_matrix(nBas,matmul(P(:,:,is),K(:,:,is))) end do - ! Total energy - - EqsGF2 = sum(ET(:)) + sum(EV(:)) + sum(EJ(:)) + sum(Ex(:)) - ! Correlation energy call UMP2(nBas,nC,nO,nV,nR,ERI_aaaa,ERI_aabb,ERI_bbbb,ENuc,EqsGF2,eGF2,Ec) + ! Total energy + + EqsGF2 = sum(ET(:)) + sum(EV(:)) + sum(EJ(:)) + sum(Ex(:)) + sum(Ec(:)) + !------------------------------------------------------------------------ ! Print results !------------------------------------------------------------------------ diff --git a/src/GF/self_energy_GF2.f90 b/src/GF/self_energy_GF2.f90 index 173ecaf..f4b8e90 100644 --- a/src/GF/self_energy_GF2.f90 +++ b/src/GF/self_energy_GF2.f90 @@ -28,7 +28,7 @@ subroutine self_energy_GF2(eta,nBas,nC,nO,nV,nR,nS,eHF,eGF2,ERI,SigC,Z) ! Initialize SigC(:,:) = 0d0 - Z(:) = 0d0 + Z(:) = 0d0 ! Compute GF2 self-energy and renormalization factor diff --git a/src/GF/unrestricted_self_energy_GF2.f90 b/src/GF/unrestricted_self_energy_GF2.f90 index 7bc83dd..cce6adf 100644 --- a/src/GF/unrestricted_self_energy_GF2.f90 +++ b/src/GF/unrestricted_self_energy_GF2.f90 @@ -36,6 +36,7 @@ subroutine unrestricted_self_energy_GF2(nBas,nC,nO,nV,nR,eta,ERI_aa,ERI_ab,ERI_b !---------------------! SigC(:,:,:) = 0d0 + Z(:,:) = 0d0 !----------------! ! Spin-up sector diff --git a/src/MP/MP2.f90 b/src/MP/MP2.f90 index 8b028df..7dd921f 100644 --- a/src/MP/MP2.f90 +++ b/src/MP/MP2.f90 @@ -23,7 +23,7 @@ subroutine MP2(nBas,nC,nO,nV,nR,ERI,ENuc,EHF,e,EcMP2) ! Output variables - double precision,intent(out) :: EcMP2(3) + double precision,intent(out) :: EcMP2 ! Hello world @@ -57,20 +57,18 @@ subroutine MP2(nBas,nC,nO,nV,nR,ERI,ENuc,EHF,e,EcMP2) enddo enddo - EcMP2(2) = 2d0*E2a - EcMP2(3) = -E2b - EcMP2(1) = EcMP2(2) + EcMP2(3) + EcMP2 = 2d0*E2a - E2b write(*,*) write(*,'(A32)') '--------------------------' write(*,'(A32)') ' MP2 calculation ' write(*,'(A32)') '--------------------------' - write(*,'(A32,1X,F16.10)') ' MP2 correlation energy = ',EcMP2(1) - write(*,'(A32,1X,F16.10)') ' Direct part = ',EcMP2(2) - write(*,'(A32,1X,F16.10)') ' Exchange part = ',EcMP2(3) + write(*,'(A32,1X,F16.10)') ' MP2 correlation energy = ',EcMP2 + write(*,'(A32,1X,F16.10)') ' Direct part = ',2d0*E2a + write(*,'(A32,1X,F16.10)') ' Exchange part = ',-E2b write(*,'(A32)') '--------------------------' - write(*,'(A32,1X,F16.10)') ' MP2 electronic energy = ',EHF + EcMP2(1) - write(*,'(A32,1X,F16.10)') ' MP2 total energy = ',ENuc + EHF + EcMP2(1) + write(*,'(A32,1X,F16.10)') ' MP2 electronic energy = ',EHF + EcMP2 + write(*,'(A32,1X,F16.10)') ' MP2 total energy = ',ENuc + EHF + EcMP2 write(*,'(A32)') '--------------------------' write(*,*)