From efaf67b6496f339209d16ffbe5e6db10de05a7a9 Mon Sep 17 00:00:00 2001 From: Pierre-Francois Loos Date: Mon, 19 Oct 2020 14:20:35 +0200 Subject: [PATCH] FC for CCD --- input/methods | 4 ++-- src/CC/CCD.f90 | 31 +++++++++++++------------- src/CC/CCD_correlation_energy.f90 | 8 +++---- src/CC/MP3_correlation_energy.f90 | 12 +++++------ src/CC/form_X.f90 | 36 +++++++++++++++---------------- src/CC/form_delta_OOVV.f90 | 10 ++++----- src/CC/form_u.f90 | 28 ++++++++++++------------ src/CC/form_v.f90 | 36 +++++++++++++++---------------- 8 files changed, 83 insertions(+), 82 deletions(-) diff --git a/input/methods b/input/methods index 2242045..d1b5eb7 100644 --- a/input/methods +++ b/input/methods @@ -1,9 +1,9 @@ # RHF UHF KS MOM T F F F # MP2* MP3 MP2-F12 - F T F + F F F # CCD CCSD CCSD(T) - F F F + T F F # drCCD rCCD lCCD pCCD F F F F # CIS* CIS(D) CID CISD diff --git a/src/CC/CCD.f90 b/src/CC/CCD.f90 index b662148..fad21b4 100644 --- a/src/CC/CCD.f90 +++ b/src/CC/CCD.f90 @@ -90,11 +90,11 @@ subroutine CCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,e ! Form energy denominator - allocate(eO(nO),eV(nV)) - allocate(delta_OOVV(nO,nO,nV,nV)) + allocate(eO(nO-nC),eV(nV-nR)) + allocate(delta_OOVV(nO-nC,nO-nC,nV-nR,nV-nR)) - eO(:) = seHF(1:nO) - eV(:) = seHF(nO+1:nBas) + eO(:) = seHF(nC+1:nO) + eV(:) = seHF(nO+1:nBas-nR) call form_delta_OOVV(nC,nO,nV,nR,eO,eV,delta_OOVV) @@ -102,18 +102,19 @@ subroutine CCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,e ! Create integral batches - allocate(OOOO(nO,nO,nO,nO),OOVV(nO,nO,nV,nV),OVOV(nO,nV,nO,nV),VVVV(nV,nV,nV,nV)) + allocate(OOOO(nO-nC,nO-nC,nO-nC,nO-nC),OOVV(nO-nC,nO-nC,nV-nR,nV-nR), & + OVOV(nO-nC,nV-nR,nO-nC,nV-nR),VVVV(nV-nR,nV-nR,nV-nR,nV-nR)) - OOOO(:,:,:,:) = dbERI( 1:nO , 1:nO , 1:nO , 1:nO ) - OOVV(:,:,:,:) = dbERI( 1:nO , 1:nO ,nO+1:nBas,nO+1:nBas) - OVOV(:,:,:,:) = dbERI( 1:nO ,nO+1:nBas, 1:nO ,nO+1:nBas) - VVVV(:,:,:,:) = dbERI(nO+1:nBas,nO+1:nBas,nO+1:nBas,nO+1:nBas) + OOOO(:,:,:,:) = dbERI(nC+1:nO ,nC+1:nO ,nC+1:nO ,nC+1:nO ) + OOVV(:,:,:,:) = dbERI(nC+1:nO ,nC+1:nO ,nO+1:nBas-nR,nO+1:nBas-nR) + OVOV(:,:,:,:) = dbERI(nC+1:nO ,nO+1:nBas-nR,nC+1:nO ,nO+1:nBas-nR) + VVVV(:,:,:,:) = dbERI(nO+1:nBas-nR,nO+1:nBas-nR,nO+1:nBas-nR,nO+1:nBas-nR) deallocate(dbERI) ! MP2 guess amplitudes - allocate(t2(nO,nO,nV,nV)) + allocate(t2(nO-nC,nO-nC,nV-nR,nV-nR)) t2(:,:,:,:) = -OOVV(:,:,:,:)/delta_OOVV(:,:,:,:) @@ -122,12 +123,12 @@ subroutine CCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,e ! Memory allocation for DIIS - allocate(error_diis(nO*nO*nV*nV,max_diis),t_diis(nO*nO*nV*nV,max_diis)) + allocate(error_diis((nO-nR)**2*(nV-nR)**2,max_diis),t_diis((nO-nR)**2*(nV-nR)**2,max_diis)) ! Initialization - allocate(r2(nO,nO,nV,nV),u(nO,nO,nV,nV),v(nO,nO,nV,nV)) - allocate(X1(nO,nO,nO,nO),X2(nV,nV),X3(nO,nO),X4(nO,nO,nV,nV)) + allocate(r2(nO-nC,nO-nC,nV-nR,nV-nR),u(nO-nC,nO-nC,nV-nR,nV-nR),v(nO-nC,nO-nC,nV-nR,nV-nR)) + allocate(X1(nO-nC,nO-nC,nO-nC,nO-nC),X2(nV-nR,nV-nR),X3(nO-nC,nO-nC),X4(nO-nC,nO-nC,nV-nR,nV-nR)) Conv = 1d0 nSCF = 0 @@ -171,7 +172,7 @@ subroutine CCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,e ! Check convergence - Conv = maxval(abs(r2(nC+1:nO,nC+1:nO,1:nV-nR,1:nV-nR))) + Conv = maxval(abs(r2(:,:,:,:))) ! Update amplitudes @@ -190,7 +191,7 @@ subroutine CCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,e ! DIIS extrapolation n_diis = min(n_diis+1,max_diis) - call DIIS_extrapolation(rcond,nO*nO*nV*nV,nO*nO*nV*nV,n_diis,error_diis,t_diis,-r2/delta_OOVV,t2) + call DIIS_extrapolation(rcond,(nO-nC)**2*(nV-nR)**2,(nO-nC)**2*(nV-nR)**2,n_diis,error_diis,t_diis,-r2/delta_OOVV,t2) ! Reset DIIS if required diff --git a/src/CC/CCD_correlation_energy.f90 b/src/CC/CCD_correlation_energy.f90 index b8bd52b..67c07a3 100644 --- a/src/CC/CCD_correlation_energy.f90 +++ b/src/CC/CCD_correlation_energy.f90 @@ -7,8 +7,8 @@ subroutine CCD_correlation_energy(nC,nO,nV,nR,OOVV,t2,EcCCD) ! Input variables integer,intent(in) :: nC,nO,nV,nR - double precision,intent(in) :: OOVV(nO,nO,nV,nV) - double precision,intent(in) :: t2(nO,nO,nV,nV) + double precision,intent(in) :: OOVV(nO-nC,nO-nC,nV-nR,nV-nR) + double precision,intent(in) :: t2(nO-nC,nO-nC,nV-nR,nV-nR) ! Local variables @@ -19,8 +19,8 @@ subroutine CCD_correlation_energy(nC,nO,nV,nR,OOVV,t2,EcCCD) double precision,intent(out) :: EcCCD EcCCD = 0d0 - do i=nC+1,nO - do j=nC+1,nO + do i=1,nO-nC + do j=1,nO-nC do a=1,nV-nR do b=1,nV-nR diff --git a/src/CC/MP3_correlation_energy.f90 b/src/CC/MP3_correlation_energy.f90 index 739335f..94f5d54 100644 --- a/src/CC/MP3_correlation_energy.f90 +++ b/src/CC/MP3_correlation_energy.f90 @@ -7,10 +7,10 @@ subroutine MP3_correlation_energy(nC,nO,nV,nR,OOVV,t2,v,delta_OOVV,EcMP3) ! Input variables integer,intent(in) :: nC,nO,nV,nR - double precision,intent(in) :: OOVV(nO,nO,nV,nV) - double precision,intent(in) :: t2(nO,nO,nV,nV) - double precision,intent(in) :: v(nO,nO,nV,nV) - double precision,intent(in) :: delta_OOVV(nO,nO,nV,nV) + double precision,intent(in) :: OOVV(nO-nC,nO-nC,nV-nR,nV-nR) + double precision,intent(in) :: t2(nO-nC,nO-nC,nV-nR,nV-nR) + double precision,intent(in) :: v(nO-nC,nO-nC,nV-nR,nV-nR) + double precision,intent(in) :: delta_OOVV(nO-nC,nO-nC,nV-nR,nV-nR) ! Local variables @@ -21,8 +21,8 @@ subroutine MP3_correlation_energy(nC,nO,nV,nR,OOVV,t2,v,delta_OOVV,EcMP3) double precision,intent(out) :: EcMP3 EcMP3 = 0d0 - do i=nC+1,nO - do j=nC+1,nO + do i=1,nO-nC + do j=1,nO-nC do a=1,nV-nR do b=1,nV-nR diff --git a/src/CC/form_X.f90 b/src/CC/form_X.f90 index 699e418..95cc7bc 100644 --- a/src/CC/form_X.f90 +++ b/src/CC/form_X.f90 @@ -7,8 +7,8 @@ subroutine form_X(nC,nO,nV,nR,OOVV,t2,X1,X2,X3,X4) ! Input variables integer,intent(in) :: nC,nO,nV,nR - double precision,intent(in) :: t2(nO,nO,nV,nV) - double precision,intent(in) :: OOVV(nO,nO,nV,nV) + double precision,intent(in) :: t2(nO-nC,nO-nC,nV-nR,nV-nR) + double precision,intent(in) :: OOVV(nO-nC,nO-nC,nV-nR,nV-nR) ! Local variables @@ -17,10 +17,10 @@ subroutine form_X(nC,nO,nV,nR,OOVV,t2,X1,X2,X3,X4) ! Output variables - double precision,intent(out) :: X1(nO,nO,nO,nO) - double precision,intent(out) :: X2(nV,nV) - double precision,intent(out) :: X3(nO,nO) - double precision,intent(out) :: X4(nO,nO,nV,nV) + double precision,intent(out) :: X1(nO-nC,nO-nC,nO-nC,nO-nC) + double precision,intent(out) :: X2(nV-nR,nV-nR) + double precision,intent(out) :: X3(nO-nC,nO-nC) + double precision,intent(out) :: X4(nO-nC,nO-nC,nV-nR,nV-nR) ! Initialization @@ -31,10 +31,10 @@ subroutine form_X(nC,nO,nV,nR,OOVV,t2,X1,X2,X3,X4) ! Build X1 - do k=nC+1,nO - do l=nC+1,nO - do i=nC+1,nO - do j=nC+1,nO + do k=1,nO-nC + do l=1,nO-nC + do i=1,nO-nC + do j=1,nO-nC do c=1,nV-nR do d=1,nV-nR X1(k,l,i,j) = X1(k,l,i,j) + OOVV(k,l,c,d)*t2(i,j,c,d) @@ -49,8 +49,8 @@ subroutine form_X(nC,nO,nV,nR,OOVV,t2,X1,X2,X3,X4) do b=1,nV-nR do c=1,nV-nR - do k=nC+1,nO - do l=nC+1,nO + do k=1,nO-nC + do l=1,nO-nC do d=1,nV-nR X2(b,c) = X2(b,c) + OOVV(k,l,c,d)*t2(k,l,b,d) enddo @@ -61,9 +61,9 @@ subroutine form_X(nC,nO,nV,nR,OOVV,t2,X1,X2,X3,X4) ! Build X3 - do k=nC+1,nO - do j=nC+1,nO - do l=nC+1,nO + do k=1,nO-nC + do j=1,nO-nC + do l=1,nO-nC do c=1,nV-nR do d=1,nV-nR X3(k,j) = X3(k,j) + OOVV(k,l,c,d)*t2(j,l,c,d) @@ -75,11 +75,11 @@ subroutine form_X(nC,nO,nV,nR,OOVV,t2,X1,X2,X3,X4) ! Build X4 - do i=nC+1,nO - do l=nC+1,nO + do i=1,nO-nC + do l=1,nO-nC do a=1,nV-nR do d=1,nV-nR - do k=nC+1,nO + do k=1,nO-nC do c=1,nV-nR X4(i,l,a,d) = X4(i,l,a,d) + OOVV(k,l,c,d)*t2(i,k,a,c) enddo diff --git a/src/CC/form_delta_OOVV.f90 b/src/CC/form_delta_OOVV.f90 index fbcd2d6..b1ecbfa 100644 --- a/src/CC/form_delta_OOVV.f90 +++ b/src/CC/form_delta_OOVV.f90 @@ -7,8 +7,8 @@ subroutine form_delta_OOVV(nC,nO,nV,nR,eO,eV,delta) ! Input variables integer,intent(in) :: nC,nO,nV,nR - double precision,intent(in) :: eO(nO) - double precision,intent(in) :: eV(nV) + double precision,intent(in) :: eO(nO-nC) + double precision,intent(in) :: eV(nV-nR) ! Local variables @@ -16,10 +16,10 @@ subroutine form_delta_OOVV(nC,nO,nV,nR,eO,eV,delta) ! Output variables - double precision,intent(out) :: delta(nO,nO,nV,nV) + double precision,intent(out) :: delta(nO-nC,nO-nC,nV-nR,nV-nR) - do i=nC+1,nO - do j=nC+1,nO + do i=1,nO-nC + do j=1,nO-nC do a=1,nV-nR do b=1,nV-nR diff --git a/src/CC/form_u.f90 b/src/CC/form_u.f90 index 739f5bd..6ae986c 100644 --- a/src/CC/form_u.f90 +++ b/src/CC/form_u.f90 @@ -7,10 +7,10 @@ subroutine form_u(nC,nO,nV,nR,OOOO,VVVV,OVOV,t2,u) ! Input variables integer,intent(in) :: nC,nO,nV,nR - double precision,intent(in) :: t2(nO,nO,nV,nV) - double precision,intent(in) :: OOOO(nO,nO,nO,nO) - double precision,intent(in) :: VVVV(nV,nV,nV,nV) - double precision,intent(in) :: OVOV(nO,nV,nO,nV) + double precision,intent(in) :: t2(nO-nC,nO-nC,nV-nR,nV-nR) + double precision,intent(in) :: OOOO(nO-nC,nO-nC,nO-nC,nO-nC) + double precision,intent(in) :: VVVV(nV-nR,nV-nR,nV-nR,nV-nR) + double precision,intent(in) :: OVOV(nO-nC,nV-nR,nO-nC,nV-nR) ! Local variables @@ -19,12 +19,12 @@ subroutine form_u(nC,nO,nV,nR,OOOO,VVVV,OVOV,t2,u) ! Output variables - double precision,intent(out) :: u(nO,nO,nV,nV) + double precision,intent(out) :: u(nO-nC,nO-nC,nV-nR,nV-nR) u(:,:,:,:) = 0d0 - do i=nC+1,nO - do j=nC+1,nO + do i=1,nO-nC + do j=1,nO-nC do a=1,nV-nR do b=1,nV-nR do c=1,nV-nR @@ -37,10 +37,10 @@ subroutine form_u(nC,nO,nV,nR,OOOO,VVVV,OVOV,t2,u) enddo enddo - do i=nC+1,nO - do j=nC+1,nO - do k=nC+1,nO - do l=nC+1,nO + do i=1,nO-nC + do j=1,nO-nC + do k=1,nO-nC + do l=1,nO-nC do a=1,nV-nR do b=1,nV-nR u(i,j,a,b) = u(i,j,a,b) + 0.5d0*OOOO(k,l,i,j)*t2(k,l,a,b) @@ -51,9 +51,9 @@ subroutine form_u(nC,nO,nV,nR,OOOO,VVVV,OVOV,t2,u) enddo enddo - do i=nC+1,nO - do j=nC+1,nO - do k=nC+1,nO + do i=1,nO-nC + do j=1,nO-nC + do k=1,nO-nC do a=1,nV-nR do b=1,nV-nR do c=1,nV-nR diff --git a/src/CC/form_v.f90 b/src/CC/form_v.f90 index 12c6150..38c839c 100644 --- a/src/CC/form_v.f90 +++ b/src/CC/form_v.f90 @@ -7,11 +7,11 @@ subroutine form_v(nC,nO,nV,nR,X1,X2,X3,X4,t2,v) ! Input variables integer,intent(in) :: nC,nO,nV,nR - double precision,intent(in) :: t2(nO,nO,nV,nV) - double precision,intent(in) :: X1(nO,nO,nO,nO) - double precision,intent(in) :: X2(nV,nV) - double precision,intent(in) :: X3(nO,nO) - double precision,intent(in) :: X4(nO,nO,nV,nV) + double precision,intent(in) :: t2(nO-nC,nO-nC,nV-nR,nV-nR) + double precision,intent(in) :: X1(nO-nC,nO-nC,nO-nC,nO-nC) + double precision,intent(in) :: X2(nV-nR,nV-nR) + double precision,intent(in) :: X3(nO-nC,nO-nC) + double precision,intent(in) :: X4(nO-nC,nO-nC,nV-nR,nV-nR) ! Local variables @@ -20,16 +20,16 @@ subroutine form_v(nC,nO,nV,nR,X1,X2,X3,X4,t2,v) ! Output variables - double precision,intent(out) :: v(nO,nO,nV,nV) + double precision,intent(out) :: v(nO-nC,nO-nC,nV-nR,nV-nR) v(:,:,:,:) = 0d0 - do i=nC+1,nO - do j=nC+1,nO + do i=1,nO-nC + do j=1,nO-nC do a=1,nV-nR do b=1,nV-nR - do k=nC+1,nO - do l=nC+1,nO + do k=1,nO-nC + do l=1,nO-nC v(i,j,a,b) = v(i,j,a,b) + 0.25d0*X1(k,l,i,j)*t2(k,l,a,b) enddo enddo @@ -38,8 +38,8 @@ subroutine form_v(nC,nO,nV,nR,X1,X2,X3,X4,t2,v) enddo enddo - do i=nC+1,nO - do j=nC+1,nO + do i=1,nO-nC + do j=1,nO-nC do a=1,nV-nR do b=1,nV-nR do c=1,nV-nR @@ -50,11 +50,11 @@ subroutine form_v(nC,nO,nV,nR,X1,X2,X3,X4,t2,v) enddo enddo - do i=nC+1,nO - do j=nC+1,nO + do i=1,nO-nC + do j=1,nO-nC do a=1,nV-nR do b=1,nV-nR - do k=nC+1,nO + do k=1,nO-nC v(i,j,a,b) = v(i,j,a,b) - 0.5d0*(X3(k,j)*t2(i,k,a,b) + X3(k,i)*t2(k,j,a,b)) enddo enddo @@ -62,11 +62,11 @@ subroutine form_v(nC,nO,nV,nR,X1,X2,X3,X4,t2,v) enddo enddo - do i=nC+1,nO - do j=nC+1,nO + do i=1,nO-nC + do j=1,nO-nC do a=1,nV-nR do b=1,nV-nR - do k=nC+1,nO + do k=1,nO-nC do c=1,nV-nR v(i,j,a,b) = v(i,j,a,b) + (X4(i,k,a,c)*t2(j,k,b,c) + X4(i,k,b,c)*t2(k,j,a,c)) enddo