4
1
mirror of https://github.com/pfloos/quack synced 2024-12-23 04:43:53 +01:00

remove construction of matrices in ppRPA

This commit is contained in:
Pierre-Francois Loos 2023-07-17 23:05:49 +02:00
parent 91954dba8d
commit 608a62dd6a
4 changed files with 84 additions and 12 deletions

View File

@ -51,6 +51,9 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,evDyn,pp
double precision :: EcAC(nspin)
double precision :: EcppBSE(nspin)
double precision :: EcGM
double precision,allocatable :: Bpp(:,:)
double precision,allocatable :: Cpp(:,:)
double precision,allocatable :: Dpp(:,:)
double precision,allocatable :: Om1ab(:),Om1aa(:)
double precision,allocatable :: X1ab(:,:),X1aa(:,:)
double precision,allocatable :: Y1ab(:,:),Y1aa(:,:)
@ -102,7 +105,15 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,evDyn,pp
! Compute linear response
call ppLR(iblock,TDA_T,nBas,nC,nO,nV,nR,nOOab,nVVab,1d0,eHF,ERI_MO,Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,EcRPA(ispin))
allocate(Bpp(nVVab,nOOab),Cpp(nVVab,nVVab),Dpp(nOOab,nOOab))
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOab,nVVab,1d0,ERI_MO,Bpp)
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nOOab,nVVab,1d0,eHF,ERI_MO,Cpp)
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOab,nVVab,1d0,eHF,ERI_MO,Dpp)
call ppLR(TDA_T,nOOab,nVVab,Bpp,Cpp,Dpp,Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,EcRPA(ispin))
deallocate(Bpp,Cpp,Dpp)
call print_excitation('pp-RPA (N+2)',iblock,nVVab,Om1ab(:))
call print_excitation('pp-RPA (N-2)',iblock,nOOab,Om2ab(:))
@ -116,7 +127,15 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,evDyn,pp
! Compute linear response
call ppLR(iblock,TDA_T,nBas,nC,nO,nV,nR,nOOaa,nVVaa,1d0,eHF,ERI_MO,Om1aa,X1aa,Y1aa,Om2aa,X2aa,Y2aa,EcRPA(ispin))
allocate(Bpp(nVVaa,nOOaa),Cpp(nVVaa,nVVaa),Dpp(nOOaa,nOOaa))
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOaa,nVVaa,1d0,ERI_MO,Bpp)
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nOOaa,nVVaa,1d0,eHF,ERI_MO,Cpp)
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOaa,nVVaa,1d0,eHF,ERI_MO,Dpp)
call ppLR(TDA_T,nOOaa,nVVaa,Bpp,Cpp,Dpp,Om1aa,X1aa,Y1aa,Om2aa,X2aa,Y2aa,EcRPA(ispin))
deallocate(Bpp,Cpp,Dpp)
call print_excitation('pp-RPA (N+2)',iblock,nVVaa,Om1aa(:))
call print_excitation('pp-RPA (N-2)',iblock,nOOaa,Om2aa(:))
@ -198,11 +217,28 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,evDyn,pp
ispin = 1
iblock = 3
call ppLR(iblock,TDA_T,nBas,nC,nO,nV,nR,nOOab,nVVab,1d0,eGT,ERI_MO,Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,EcRPA(ispin))
allocate(Bpp(nVVab,nOOab),Cpp(nVVab,nVVab),Dpp(nOOab,nOOab))
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOab,nVVab,1d0,ERI_MO,Bpp)
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nOOab,nVVab,1d0,eGT,ERI_MO,Cpp)
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOab,nVVab,1d0,eGT,ERI_MO,Dpp)
call ppLR(TDA_T,nOOab,nVVab,Bpp,Cpp,Dpp,Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,EcRPA(ispin))
deallocate(Bpp,Cpp,Dpp)
ispin = 2
iblock = 4
call ppLR(iblock,TDA_T,nBas,nC,nO,nV,nR,nOOaa,nVVaa,1d0,eGT,ERI_MO,Om1aa,X1aa,Y1aa,Om2aa,X2aa,Y2aa,EcRPA(ispin))
allocate(Bpp(nVVaa,nOOaa),Cpp(nVVaa,nVVaa),Dpp(nOOaa,nOOaa))
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOaa,nVVaa,1d0,ERI_MO,Bpp)
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nOOaa,nVVaa,1d0,eGT,ERI_MO,Cpp)
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOaa,nVVaa,1d0,eGT,ERI_MO,Dpp)
call ppLR(TDA_T,nOOaa,nVVaa,Bpp,Cpp,Dpp,Om1aa,X1aa,Y1aa,Om2aa,X2aa,Y2aa,EcRPA(ispin))
deallocate(Bpp,Cpp,Dpp)
EcRPA(1) = EcRPA(1) - EcRPA(2)
EcRPA(2) = 3d0*EcRPA(2)

View File

@ -60,6 +60,9 @@ subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T
double precision,allocatable :: e_diis(:,:)
double precision,allocatable :: eGT(:)
double precision,allocatable :: eOld(:)
double precision,allocatable :: Bpp(:,:)
double precision,allocatable :: Cpp(:,:)
double precision,allocatable :: Dpp(:,:)
double precision,allocatable :: Om1s(:),Om1t(:)
double precision,allocatable :: X1s(:,:),X1t(:,:)
double precision,allocatable :: Y1s(:,:),Y1t(:,:)
@ -132,7 +135,15 @@ subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T
! Compute linear response
call ppLR(iblock,TDA_T,nBas,nC,nO,nV,nR,nOOs,nVVs,1d0,eGT,ERI_MO,Om1s,X1s,Y1s,Om2s,X2s,Y2s,EcRPA(ispin))
allocate(Bpp(nVVs,nOOs),Cpp(nVVs,nVVs),Dpp(nOOs,nOOs))
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,1d0,ERI_MO,Bpp)
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,1d0,eHF,ERI_MO,Cpp)
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,1d0,eHF,ERI_MO,Dpp)
call ppLR(TDA_T,nOOs,nVVs,Bpp,Cpp,Dpp,Om1s,X1s,Y1s,Om2s,X2s,Y2s,EcRPA(ispin))
deallocate(Bpp,Cpp,Dpp)
!----------------------------------------------
! alpha-alpha block
@ -143,7 +154,15 @@ subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T
! Compute linear response
call ppLR(iblock,TDA_T,nBas,nC,nO,nV,nR,nOOt,nVVt,1d0,eGT,ERI_MO,Om1t,X1t,Y1t,Om2t,X2t,Y2t,EcRPA(ispin))
allocate(Bpp(nVVt,nOOt),Cpp(nVVt,nVVt),Dpp(nOOt,nOOt))
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,1d0,ERI_MO,Bpp)
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,1d0,eHF,ERI_MO,Cpp)
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,1d0,eHF,ERI_MO,Dpp)
call ppLR(TDA_T,nOOt,nVVt,Bpp,Cpp,Dpp,Om1t,X1t,Y1t,Om2t,X2t,Y2t,EcRPA(ispin))
deallocate(Bpp,Cpp,Dpp)
EcRPA(1) = EcRPA(1) - EcRPA(2)
EcRPA(2) = 3d0*EcRPA(2)

View File

@ -76,6 +76,9 @@ subroutine qsGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T
double precision,allocatable :: cp(:,:)
double precision,allocatable :: eGT(:)
double precision,allocatable :: eOld(:)
double precision,allocatable :: Bpp(:,:)
double precision,allocatable :: Cpp(:,:)
double precision,allocatable :: Dpp(:,:)
double precision,allocatable :: Om1s(:),Om1t(:)
double precision,allocatable :: X1s(:,:),X1t(:,:)
double precision,allocatable :: Y1s(:,:),Y1t(:,:)
@ -188,12 +191,28 @@ subroutine qsGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T
ispin = 1
iblock = 3
call ppLR(iblock,TDA_T,nBas,nC,nO,nV,nR,nOOs,nVVs,1d0,eGT,ERI_MO,Om1s,X1s,Y1s,Om2s,X2s,Y2s,EcRPA(ispin))
allocate(Bpp(nVVs,nOOs),Cpp(nVVs,nVVs),Dpp(nOOs,nOOs))
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,1d0,ERI_MO,Bpp)
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,1d0,eGT,ERI_MO,Cpp)
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,1d0,eGT,ERI_MO,Dpp)
call ppLR(TDA_T,nOOs,nVVs,Bpp,Cpp,Dpp,Om1s,X1s,Y1s,Om2s,X2s,Y2s,EcRPA(ispin))
deallocate(Bpp,Cpp,Dpp)
ispin = 2
iblock = 4
call ppLR(iblock,TDA_T,nBas,nC,nO,nV,nR,nOOt,nVVt,1d0,eGT,ERI_MO,Om1t,X1t,Y1t,Om2t,X2t,Y2t,EcRPA(ispin))
allocate(Bpp(nVVt,nOOt),Cpp(nVVt,nVVt),Dpp(nOOt,nOOt))
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,1d0,ERI_MO,Bpp)
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,1d0,eGT,ERI_MO,Cpp)
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,1d0,eGT,ERI_MO,Dpp)
call ppLR(TDA_T,nOOt,nVVt,Bpp,Cpp,Dpp,Om1t,X1t,Y1t,Om2t,X2t,Y2t,EcRPA(ispin))
deallocate(Bpp,Cpp,Dpp)
EcRPA(1) = EcRPA(1) - EcRPA(2)
EcRPA(2) = 3d0*EcRPA(2)

View File

@ -1,4 +1,4 @@
subroutine ppLR(TDA,nBas,nOO,nVV,Om1,X1,Y1,Om2,X2,Y2,EcRPA)
subroutine ppLR(TDA,nOO,nVV,B,C,D,Om1,X1,Y1,Om2,X2,Y2,EcRPA)
! Compute the p-p channel of the linear response: see Scuseria et al. JCP 139, 104113 (2013)
@ -61,8 +61,6 @@ subroutine ppLR(TDA,nBas,nOO,nVV,Om1,X1,Y1,Om2,X2,Y2,EcRPA)
else
call ppLR_B(ispin,nBas,nC,nO,nV,nR,nOO,nVV,lambda,ERI,B)
! Diagonal blocks
M( 1:nVV , 1:nVV) = + C(1:nVV,1:nVV)