PTEROSOR/QuAcK
10
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2024-12-22 20:34:46 +01:00
Code Releases Activity

DIP- and DEA-EOM-CCD

Browse Source
This commit is contained in:
Pierre-Francois Loos 2022-09-27 14:07:31 +02:00
parent 68008e5738
commit 1f2f6634bb
11 changed files with 333 additions and 54 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
input/methods
View File

@ -3,7 +3,7 @@
# MP2* MP3 MP2-F12
F F F
# CCD pCCD DCD CCSD CCSD(T)
F F F F F
T F F F F
# drCCD rCCD crCCD lCCD
F F F F
# CIS* CIS(D) CID CISD FCI
@ -15,7 +15,7 @@
# G0W0* evGW* qsGW* ufG0W0 ufGW
F F F F F
# G0T0 evGT qsGT
T F F
F F F
# MCMP2
F
# * unrestricted version available

4
input/options
View File

@ -5,7 +5,7 @@
# CC: maxSCF thresh DIIS n_diis
64 0.00001 T 5
# spin: TDA singlet triplet spin_conserved spin_flip
T T F T T
T T T T T
# GF: maxSCF thresh DIIS n_diis lin eta renorm reg
256 0.00001 T 5 T 0.0 3 F
# GW: maxSCF thresh DIIS n_diis lin eta COHSEX SOSEX TDA_W G0W GW0 reg
@ -15,6 +15,6 @@
# ACFDT: AC Kx XBS
F F T
# BSE: BSE dBSE dTDA evDyn ppBSE
F F T F T
T F T F F
# MCMP2: nMC nEq nWalk dt nPrint iSeed doDrift
1000000 100000 10 0.3 10000 1234 T

60
src/CC/CCD.f90
View File

@ -52,14 +52,20 @@ subroutine CCD(BSE,maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ER
double precision,allocatable :: u(:,:,:,:)
double precision,allocatable :: v(:,:,:,:)
double precision,allocatable :: r2(:,:,:,:)
double precision,allocatable :: t2(:,:,:,:)
double precision,allocatable :: r(:,:,:,:)
double precision,allocatable :: t(:,:,:,:)
integer :: n_diis,i,j,a,b
double precision :: rcond
double precision,allocatable :: error_diis(:,:)
double precision,allocatable :: t_diis(:,:)
logical :: do_EE_EOM_CC_1h1p = .false.
logical :: do_EA_EOM_CC_1p = .false.
logical :: do_IP_EOM_CC_1h = .false.
logical :: do_DEA_EOM_CC_2p = .true.
logical :: do_DIP_EOM_CC_2h = .true.
! Hello world
write(*,*)
@ -123,11 +129,11 @@ subroutine CCD(BSE,maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ER
! MP2 guess amplitudes
allocate(t2(nO-nC,nO-nC,nV-nR,nV-nR))
allocate(t(nO-nC,nO-nC,nV-nR,nV-nR))
t2(:,:,:,:) = -OOVV(:,:,:,:)/delta_OOVV(:,:,:,:)
t(:,:,:,:) = -OOVV(:,:,:,:)/delta_OOVV(:,:,:,:)
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t2,EcMP2)
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t,EcMP2)
EcMP4 = 0d0
! Memory allocation for DIIS
@ -136,7 +142,7 @@ subroutine CCD(BSE,maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ER
! Initialization
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(r(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
@ -165,33 +171,33 @@ subroutine CCD(BSE,maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ER
! Form linear array
call form_u(nC,nO,nV,nR,OOOO,VVVV,OVOV,t2,u)
call form_u(nC,nO,nV,nR,OOOO,VVVV,OVOV,t,u)
! Form interemediate arrays
call form_X(nC,nO,nV,nR,OOVV,t2,X1,X2,X3,X4)
call form_X(nC,nO,nV,nR,OOVV,t,X1,X2,X3,X4)
! Form quadratic array
call form_v(nC,nO,nV,nR,X1,X2,X3,X4,t2,v)
call form_v(nC,nO,nV,nR,X1,X2,X3,X4,t,v)
! Compute residual
r2(:,:,:,:) = OOVV(:,:,:,:) + delta_OOVV(:,:,:,:)*t2(:,:,:,:) + u(:,:,:,:) + v(:,:,:,:)
r(:,:,:,:) = OOVV(:,:,:,:) + delta_OOVV(:,:,:,:)*t(:,:,:,:) + u(:,:,:,:) + v(:,:,:,:)
! Check convergence
Conv = maxval(abs(r2(:,:,:,:)))
Conv = maxval(abs(r(:,:,:,:)))
! Update amplitudes
t2(:,:,:,:) = t2(:,:,:,:) - r2(:,:,:,:)/delta_OOVV(:,:,:,:)
t(:,:,:,:) = t(:,:,:,:) - r(:,:,:,:)/delta_OOVV(:,:,:,:)
! Compute correlation energy
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t2,EcCCD)
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t,EcCCD)
if(nSCF == 1) call MP3_correlation_energy(nC,nO,nV,nR,OOVV,t2,v,delta_OOVV,EcMP3)
if(nSCF == 1) call MP3_correlation_energy(nC,nO,nV,nR,OOVV,t,v,delta_OOVV,EcMP3)
! Dump results
@ -200,7 +206,7 @@ subroutine CCD(BSE,maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ER
! DIIS extrapolation
n_diis = min(n_diis+1,max_diis)
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)
call DIIS_extrapolation(rcond,(nO-nC)**2*(nV-nR)**2,(nO-nC)**2*(nV-nR)**2,n_diis,error_diis,t_diis,-r/delta_OOVV,t)
! Reset DIIS if required
@ -246,4 +252,28 @@ subroutine CCD(BSE,maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ER
write(*,'(1X,A15,1X,F10.6)') 'Ec(MP4-SDQ) = ',EcMP4
write(*,*)
!------------------------------------------------------------------------
! EOM section
!------------------------------------------------------------------------
! EE-EOM-CCD (1h1p)
! if(do_EE-EOM-CC_1h1p) call EE-EOM-CCD_1h1p()
! EA-EOM (1p)
! if(do_EA-EOM-CC_1p) call EA-EOM-CCD_1p()
! IP-EOM-CCD(1h)
! if(do_IP-EOM-CC_1h) call IP-EOM-CCD_1h()
! DEA-EOM (2p)
if(do_DEA_EOM_CC_2p) call DEA_EOM_CCD_2p(nC,nO,nV,nR,eV,OOVV,VVVV,t)
! DIP-EOM-CCD(2h)
if(do_DIP_EOM_CC_2h) call DIP_EOM_CCD_2h(nC,nO,nV,nR,eO,OOVV,OOOO,t)
end subroutine CCD

116
src/CC/DEA_EOM_CCD_2p.f90 Normal file
View File

@ -0,0 +1,116 @@
subroutine DEA_EOM_CCD_2p(nC,nO,nV,nR,eV,OOVV,VVVV,t)
! DEA-EOM-CCD calculation up to 2p
implicit none
! Input variables
integer,intent(in) :: nC
integer,intent(in) :: nO
integer,intent(in) :: nV
integer,intent(in) :: nR
double precision,intent(in) :: eV(nV)
double precision,intent(in) :: OOVV(nO,nO,nV,nV)
double precision,intent(in) :: VVVV(nV,nV,nV,nV)
double precision,intent(in) :: t(nO,nO,nV,nV)
! Local variables
integer :: a,b,c,d,ab,cd
integer :: i,j
integer :: nVV
double precision,external :: Kronecker_delta
double precision,allocatable :: F(:,:)
double precision,allocatable :: W(:,:,:,:)
double precision,allocatable :: H(:,:)
double precision,allocatable :: Om(:)
! Hello world
write(*,*)
write(*,*)'********************'
write(*,*)'| DEA-EOM-CCD (2p) |'
write(*,*)'********************'
write(*,*)
! Size of the EOM Hamiltonian
nVV = nV*(nV-1)/2
! Memory allocation
allocate(F(nV,nV),W(nV,nV,nV,nV),H(nVV,nVV),Om(nVV))
! Form one-body terms
do a=1,nV-nR
do b=1,nV-nR
F(a,b) = eV(a)*Kronecker_delta(a,b)
do i=1,nO-nR
do j=1,nO-nR
do c=1,nV-nC
F(a,b) = F(a,b) - 0.5d0*OOVV(i,j,b,c)*t(i,j,a,c)
end do
end do
end do
end do
end do
! Form two-body terms
do a=1,nV-nR
do b=1,nV-nR
do c=1,nV-nR
do d=1,nV-nR
W(a,b,c,d) = VVVV(a,b,c,d)
do i=1,nO-nC
do j=i+1,nO-nC
W(a,b,c,d) = W(a,b,c,d) + OOVV(i,j,a,b)*t(i,j,a,b)
end do
end do
end do
end do
end do
end do
! Form EOM Hamiltonian
ab = 0
do a=1,nV-nR
do b=a+1,nV-nR
ab = ab + 1
cd = 0
do c=1,nV-nR
do d=c+1,nV-nR
cd = cd + 1
H(ab,cd) = F(a,c)*Kronecker_delta(b,d) + Kronecker_delta(a,c)*F(b,d) + W(a,b,c,d)
end do
end do
end do
end do
! Diagonalize EOM Hamiltonian
if(nVV > 0) call diagonalize_matrix(nVV,H,Om)
! Dump results
call print_excitation('DEA-EOM-CCD ',3,nVV,Om)
end subroutine DEA_EOM_CCD_2p

116
src/CC/DIP_EOM_CCD_2h.f90 Normal file
View File

@ -0,0 +1,116 @@
subroutine DIP_EOM_CCD_2h(nC,nO,nV,nR,eO,OOVV,OOOO,t)
! DIP-EOM-CCD calculation up to 2h
implicit none
! Input variables
integer,intent(in) :: nC
integer,intent(in) :: nO
integer,intent(in) :: nV
integer,intent(in) :: nR
double precision,intent(in) :: eO(nO)
double precision,intent(in) :: OOVV(nO,nO,nV,nV)
double precision,intent(in) :: OOOO(nO,nO,nO,nO)
double precision,intent(in) :: t(nO,nO,nV,nV)
! Local variables
integer :: i,j,k,l,ij,kl
integer :: a,b
integer :: nOO
double precision,external :: Kronecker_delta
double precision,allocatable :: F(:,:)
double precision,allocatable :: W(:,:,:,:)
double precision,allocatable :: H(:,:)
double precision,allocatable :: Om(:)
! Hello world
write(*,*)
write(*,*)'********************'
write(*,*)'| DIP-EOM-CCD (2h) |'
write(*,*)'********************'
write(*,*)
! Size of the EOM Hamiltonian
nOO = nO*(nO-1)/2
! Memory allocation
allocate(F(nO,nO),W(nO,nO,nO,nO),H(nOO,nOO),Om(nOO))
! Form one-body terms
do i=1,nO-nR
do j=1,nO-nR
F(i,j) = eO(i)*Kronecker_delta(i,j)
do k=1,nO-nR
do a=1,nV-nC
do b=1,nV-nC
F(i,j) = F(i,j) + 0.5d0*OOVV(i,k,a,b)*t(j,k,a,b)
end do
end do
end do
end do
end do
! Form two-body terms
do i=1,nO-nC
do j=1,nO-nC
do k=1,nO-nC
do l=1,nO-nC
W(i,j,k,l) = OOOO(i,j,k,l)
do a=1,nV-nR
do b=a+1,nV-nR
W(i,j,k,l) = W(i,j,k,l) + OOVV(i,j,a,b)*t(i,j,a,b)
end do
end do
end do
end do
end do
end do
! Form EOM Hamiltonian
ij = 0
do i=1,nO-nC
do j=i+1,nO-nC
ij = ij + 1
kl = 0
do k=1,nO-nC
do l=k+1,nO-nC
kl = kl + 1
H(ij,kl) = - F(i,k)*Kronecker_delta(j,l) - Kronecker_delta(i,k)*F(j,l) + W(i,j,k,l)
end do
end do
end do
end do
! Diagonalize EOM Hamiltonian
if(nOO > 0) call diagonalize_matrix(nOO,H,Om)
! Dump results
call print_excitation('DIP-EOM-CCD ',3,nOO,Om)
end subroutine DIP_EOM_CCD_2h

13
src/CC/rCCD.f90
View File

@ -1,4 +1,4 @@
subroutine rCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
subroutine rCCD(BSE,maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
! Ring CCD module
@ -6,6 +6,7 @@ subroutine rCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
! Input variables
logical,intent(in) :: BSE
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
@ -74,7 +75,15 @@ subroutine rCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
allocate(dbERI(nBas,nBas,nBas,nBas))
call antisymmetrize_ERI(2,nBas,sERI,dbERI)
if(BSE) then
call static_screening(nBas,nC,nO,nV,nR,seHF,sERI,dbERI)
else
call antisymmetrize_ERI(2,nBas,sERI,dbERI)
end if
deallocate(sERI)

8
src/GT/static_Tmatrix_C_pp.f90
View File

@ -55,15 +55,17 @@ subroutine static_Tmatrix_C_pp(ispin,eta,nBas,nC,nO,nV,nR,nOO,nVV,nOOx,nVVx,lamb
do ef=1,nVV
eps = + Om1(ef)
chi = chi + rho1(a,b,ef)*rho1(c,d,ef)*eps/(eps**2 + eta**2)
chi = chi + rho1(a,b,ef)*rho1(c,d,ef)*eps/(eps**2 + eta**2) &
+ rho1(a,b,ef)*rho1(d,c,ef)*eps/(eps**2 + eta**2)
end do
do mn=1,nOO
eps = - Om2(mn)
chi = chi + rho2(a,b,mn)*rho2(c,d,mn)*eps/(eps**2 + eta**2)
chi = chi + rho2(a,b,mn)*rho2(c,d,mn)*eps/(eps**2 + eta**2) &
+ rho2(a,b,mn)*rho2(d,c,mn)*eps/(eps**2 + eta**2)
end do
TC(ab,cd) = lambda*chi/sqrt((1d0 + Kronecker_delta(a,b))*(1d0 + Kronecker_delta(c,d)))
TC(ab,cd) = 0.5d0*lambda*chi/sqrt((1d0 + Kronecker_delta(a,b))*(1d0 + Kronecker_delta(c,d)))
end do
end do

2
src/GW/static_screening_WC_pp.f90
View File

@ -56,7 +56,7 @@ subroutine static_screening_WC_pp(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,E
do m=1,nS
eps = Omega(m)**2 + eta**2
chi = chi + rho(a,c,m)*rho(b,d,m)*Omega(m)/eps &
- rho(a,d,m)*rho(b,c,m)*Omega(m)/eps
+ rho(a,d,m)*rho(b,c,m)*Omega(m)/eps
enddo
WC(ab,cd) = + 4d0*lambda*chi/sqrt((1d0 + Kronecker_delta(a,b))*(1d0 + Kronecker_delta(c,d)))

2
src/LR/print_excitation.f90
View File

@ -33,7 +33,7 @@ subroutine print_excitation(method,ispin,nS,Omega)
'|','State','|',' Excitation energy (au) ','|',' Excitation energy (eV) ','|'
write(*,*)'-------------------------------------------------------------'
do ia=1,maxS
do ia=1,min(maxS,nS)
write(*,'(1X,A1,1X,I5,1X,A1,1X,F23.6,1X,A1,1X,F23.6,1X,A1,1X)') &
'|',ia,'|',Omega(ia),'|',Omega(ia)*HaToeV,'|'
enddo

50
src/QuAcK/QuAcK.f90
View File

@ -582,15 +582,15 @@ program QuAcK
write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for CCD = ',t_CCD,' seconds'
write(*,*)
call cpu_time(start_CCD)
call G0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,evDyn,ppBSE,singlet,triplet, &
linGW,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,eHF,Vxc,eG0W0)
call CCD(.true.,maxSCF_CC,thresh_CC,n_diis_CC,nBas,nC,nO,nV,nR,ERI_MO,ENuc,ERHF,eG0W0)
call cpu_time(end_CCD)
! call cpu_time(start_CCD)
! call G0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,evDyn,ppBSE,singlet,triplet, &
! linGW,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,eHF,Vxc,eG0W0)
! call CCD(.true.,maxSCF_CC,thresh_CC,n_diis_CC,nBas,nC,nO,nV,nR,ERI_MO,ENuc,ERHF,eG0W0)
! call cpu_time(end_CCD)
t_CCD = end_CCD - start_CCD
write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for CCD@BSE = ',t_CCD,' seconds'
write(*,*)
! t_CCD = end_CCD - start_CCD
! write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for CCD@BSE = ',t_CCD,' seconds'
! write(*,*)
end if
@ -627,16 +627,16 @@ program QuAcK
write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for CCSD or CCSD(T)= ',t_CCSD,' seconds'
write(*,*)
call cpu_time(start_CCSD)
call G0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,evDyn,ppBSE,singlet,triplet, &
linGW,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,eHF,Vxc,eG0W0)
! call cpu_time(start_CCSD)
! call G0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,evDyn,ppBSE,singlet,triplet, &
! linGW,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,eHF,Vxc,eG0W0)
call CCSD(.true.,maxSCF_CC,thresh_CC,n_diis_CC,doCCSDT,nBas,nC,nO,nV,nR,ERI_MO,ENuc,ERHF,eG0W0)
call cpu_time(end_CCSD)
! call CCSD(.true.,maxSCF_CC,thresh_CC,n_diis_CC,doCCSDT,nBas,nC,nO,nV,nR,ERI_MO,ENuc,ERHF,eG0W0)
! call cpu_time(end_CCSD)
t_CCSD = end_CCSD - start_CCSD
write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for CCSD or CCSD(T)= ',t_CCSD,' seconds'
write(*,*)
! t_CCSD = end_CCSD - start_CCSD
! write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for CCSD or CCSD(T)= ',t_CCSD,' seconds'
! write(*,*)
end if
@ -647,8 +647,7 @@ program QuAcK
if(do_drCCD) then
call cpu_time(start_CCD)
call drCCD(maxSCF_CC,thresh_CC,n_diis_CC,nBas,nC,nO,nV,nR, &
ERI_MO,ENuc,ERHF,eHF)
call drCCD(maxSCF_CC,thresh_CC,n_diis_CC,nBas,nC,nO,nV,nR,ERI_MO,ENuc,ERHF,eHF)
call cpu_time(end_CCD)
t_CCD = end_CCD - start_CCD
@ -664,8 +663,16 @@ program QuAcK
if(do_rCCD) then
call cpu_time(start_CCD)
call rCCD(maxSCF_CC,thresh_CC,n_diis_CC,nBas,nC,nO,nV,nR, &
ERI_MO,ENuc,ERHF,eHF)
call rCCD(.false.,maxSCF_CC,thresh_CC,n_diis_CC,nBas,nC,nO,nV,nR,ERI_MO,ENuc,ERHF,eHF)
call cpu_time(end_CCD)
t_CCD = end_CCD - start_CCD
write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for ring CCD = ',t_CCD,' seconds'
write(*,*)
call cpu_time(start_CCD)
call G0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,evDyn,ppBSE,singlet,triplet, &
linGW,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,eHF,Vxc,eG0W0)
call rCCD(.true.,maxSCF_CC,thresh_CC,n_diis_CC,nBas,nC,nO,nV,nR,ERI_MO,ENuc,ERHF,eG0W0)
call cpu_time(end_CCD)
t_CCD = end_CCD - start_CCD
@ -681,8 +688,7 @@ program QuAcK
if(do_crCCD) then
call cpu_time(start_CCD)
call crCCD(maxSCF_CC,thresh_CC,n_diis_CC,nBas,nC,nO,nV,nR, &
ERI_MO,ENuc,ERHF,eHF)
call crCCD(maxSCF_CC,thresh_CC,n_diis_CC,nBas,nC,nO,nV,nR,ERI_MO,ENuc,ERHF,eHF)
call cpu_time(end_CCD)

12
src/RPA/ppRPA.f90
View File

@ -74,10 +74,10 @@ subroutine ppRPA(TDA,doACFDT,singlet,triplet,nBas,nC,nO,nV,nR,ENuc,ERHF,ERI,dipo
call linear_response_pp(ispin,TDA,nBas,nC,nO,nV,nR,nOO,nVV,1d0,e,ERI, &
Omega1,X1,Y1,Omega2,X2,Y2,Ec_ppRPA(ispin))
call print_transition_vectors_pp(.true.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Omega1,X1,Y1,Omega2,X2,Y2)
! call print_transition_vectors_pp(.true.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Omega1,X1,Y1,Omega2,X2,Y2)
! call print_excitation('pp-BSE (N+2)',ispin,nVV,Omega1)
! call print_excitation('pp-BSE (N-2)',ispin,nOO,Omega2)
call print_excitation('pp-BSE (N+2)',ispin,nVV,Omega1)
call print_excitation('pp-BSE (N-2)',ispin,nOO,Omega2)
deallocate(Omega1,X1,Y1,Omega2,X2,Y2)
@ -102,10 +102,10 @@ subroutine ppRPA(TDA,doACFDT,singlet,triplet,nBas,nC,nO,nV,nR,ENuc,ERHF,ERI,dipo
call linear_response_pp(ispin,TDA,nBas,nC,nO,nV,nR,nOO,nVV,1d0,e,ERI, &
Omega1,X1,Y1,Omega2,X2,Y2,Ec_ppRPA(ispin))
call print_transition_vectors_pp(.false.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Omega1,X1,Y1,Omega2,X2,Y2)
! call print_transition_vectors_pp(.false.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Omega1,X1,Y1,Omega2,X2,Y2)
! call print_excitation('pp-BSE (N+2)',ispin,nVV,Omega1)
! call print_excitation('pp-BSE (N-2)',ispin,nOO,Omega2)
call print_excitation('pp-BSE (N+2)',ispin,nVV,Omega1)
call print_excitation('pp-BSE (N-2)',ispin,nOO,Omega2)
deallocate(Omega1,X1,Y1,Omega2,X2,Y2)