LCPQ/quack
4
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2024-12-22 20:35:36 +01:00
Code Releases Activity

GMP2

Browse Source
This commit is contained in:
Pierre-Francois Loos 2023-10-26 22:24:22 +02:00
parent 5bf12265b5
commit 7e08e6abbd
3 changed files with 27 additions and 27 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
input/methods
View File

@ -1,7 +1,7 @@
# RHF UHF GHF ROHF
F F T F
# MP2 MP3
F F
T F
# CCD pCCD DCD CCSD CCSD(T)
F F F F F
# drCCD rCCD crCCD lCCD

34
src/MP/GMP2.f90
View File

@ -22,6 +22,7 @@ subroutine GMP2(regularize,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,e)
integer :: i,j,a,b
double precision :: kappa,sigm1,sigm2
double precision :: num
double precision :: Dijab
double precision :: fs,fs2,fk
@ -76,28 +77,29 @@ subroutine GMP2(regularize,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,e)
fs2 = (1d0 - exp(-sigm2*Dijab*Dijab))/Dijab
fk = (1d0 - exp(-kappa*Dijab))**2/Dijab
E2d = E2d - 0.25d0*(ERI(i,j,a,b)-ERI(i,j,b,a))**2/Dijab
E2ds = E2ds - ERI(i,j,a,b)*ERI(i,j,a,b)*fs
E2ds2 = E2ds2 - ERI(i,j,a,b)*ERI(i,j,a,b)*fs2
E2dk = E2dk - ERI(i,j,a,b)*ERI(i,j,a,b)*fk
num = 0.25d0*(ERI(i,j,a,b)**2 + ERI(i,j,b,a)**2)
E2d = E2d - num/Dijab
E2ds = E2ds - num*fs
E2ds2 = E2ds2 - num*fs2
E2dk = E2dk - num*fk
! Second-order exchange diagram
E2x = E2x - ERI(i,j,a,b)*ERI(i,j,b,a)/Dijab
E2x = 0d0
E2xs = E2xs - ERI(i,j,a,b)*ERI(i,j,b,a)*fs
E2xs2 = E2xs2 - ERI(i,j,a,b)*ERI(i,j,b,a)*fs2
E2xk = E2xk - ERI(i,j,a,b)*ERI(i,j,b,a)*fk
num = 0.5d0*ERI(i,j,a,b)*ERI(i,j,b,a)
E2x = E2x + num/Dijab
E2xs = E2xs + num*fs
E2xs2 = E2xs2 + num*fs2
E2xk = E2xk + num*fk
enddo
enddo
enddo
enddo
EcMP2 = E2d - E2x
EcsMP2 = E2ds - E2xs
Ecs2MP2 = E2ds2 - E2xs2
EckMP2 = E2dk - E2xk
EcMP2 = E2d + E2x
EcsMP2 = E2ds + E2xs
Ecs2MP2 = E2ds2 + E2xs2
EckMP2 = E2dk + E2xk
!------------!
! MP2 energy !
@ -109,7 +111,7 @@ subroutine GMP2(regularize,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,e)
write(*,'(A32)') '---------------------------'
write(*,'(A32,1X,F16.10)') ' GMP2 correlation energy = ',EcMP2
write(*,'(A32,1X,F16.10)') ' Direct part = ',E2d
write(*,'(A32,1X,F16.10)') ' Exchange part = ',-E2x
write(*,'(A32,1X,F16.10)') ' Exchange part = ',E2x
write(*,'(A32)') '---------------------------'
write(*,'(A32,1X,F16.10)') ' GMP2 electronic energy = ',EHF + EcMP2
write(*,'(A32,1X,F16.10)') ' GMP2 total energy = ',ENuc + EHF + EcMP2
@ -128,7 +130,7 @@ subroutine GMP2(regularize,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,e)
write(*,'(A32)') '---------------------------'
write(*,'(A32,1X,F16.10)') ' GMP2 correlation energy = ',EcsMP2
write(*,'(A32,1X,F16.10)') ' Direct part = ',E2ds
write(*,'(A32,1X,F16.10)') ' Exchange part = ',-E2xs
write(*,'(A32,1X,F16.10)') ' Exchange part = ',E2xs
write(*,'(A32)') '---------------------------'
write(*,'(A32,1X,F16.10)') ' GMP2 electronic energy = ',EHF + EcsMP2
write(*,'(A32,1X,F16.10)') ' GMP2 total energy = ',ENuc + EHF + EcsMP2
@ -145,7 +147,7 @@ subroutine GMP2(regularize,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,e)
write(*,'(A32)') '---------------------------'
write(*,'(A32,1X,F16.10)') ' GMP2 correlation energy = ',Ecs2MP2
write(*,'(A32,1X,F16.10)') ' Direct part = ',E2ds2
write(*,'(A32,1X,F16.10)') ' Exchange part = ',-E2xs2
write(*,'(A32,1X,F16.10)') ' Exchange part = ',E2xs2
write(*,'(A32)') '---------------------------'
write(*,'(A32,1X,F16.10)') ' GMP2 electronic energy = ',EHF + Ecs2MP2
write(*,'(A32,1X,F16.10)') ' GMP2 total energy = ',ENuc + EHF + Ecs2MP2

18
src/QuAcK/GQuAcK.f90
View File

@ -120,21 +120,19 @@ subroutine GQuAcK(doGHF,dostab,doMP2,dophRPA,dophRPAx,doppRPA,doG0W0,doevGW,doqs
write(*,*) 'AO to MO transformation... Please be patient'
write(*,*)
! Read and transform dipole-related integrals
dipole_int_MO(:,:,:) = dipole_int_AO(:,:,:)
do ixyz=1,ncart
call AOtoMO_transform(nBas2,cHF,dipole_int_MO(:,:,ixyz))
end do
! 4-index transform
allocate(Ca(nBas,nBas2),Cb(nBas,nBas2),ERI_tmp(nBas2,nBas2,nBas2,nBas2))
Ca(:,:) = cHF(1:nBas,1:nBas2)
Cb(:,:) = cHF(nBas+1:nBas2,1:nBas2)
! Transform dipole-related integrals
do ixyz=1,ncart
call AOtoMO_transform_GHF(nBas,nBas2,Ca,Cb,dipole_int_AO(:,:,ixyz),dipole_int_MO(:,:,ixyz))
end do
! 4-index transform
call AOtoMO_integral_transform_GHF(nBas,nBas2,Ca,Ca,Ca,Ca,ERI_AO,ERI_tmp)
ERI_MO(:,:,:,:) = ERI_tmp(:,:,:,:)