LCPQ/quack
4
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2024-11-20 13:02:38 +01:00
Code Releases Activity

rename LR routines

Browse Source
This commit is contained in:
Pierre-Francois Loos 2023-07-28 14:14:35 +02:00
parent 9c7725a0d9
commit 4576ee06c3
35 changed files with 129 additions and 366 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

236
src/CC/BCCD.f90
View File

@ -1,236 +0,0 @@
subroutine BCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
! Brueckner CCD module
implicit none
! Input variables
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
integer,intent(in) :: nBasin
integer,intent(in) :: nCin
integer,intent(in) :: nOin
integer,intent(in) :: nVin
integer,intent(in) :: nRin
double precision,intent(in) :: ENuc,ERHF
double precision,intent(in) :: eHF(nBasin)
double precision,intent(in) :: ERI(nBasin,nBasin,nBasin,nBasin)
! Local variables
integer :: nBas
integer :: nC
integer :: nO
integer :: nV
integer :: nR
integer :: nSCF
double precision :: Conv
double precision :: EcMP2,EcMP3,EcMP4
double precision :: ECCD,EcCCD
double precision,allocatable :: seHF(:)
double precision,allocatable :: sERI(:,:,:,:)
double precision,allocatable :: dbERI(:,:,:,:)
double precision,allocatable :: eO(:)
double precision,allocatable :: eV(:)
double precision,allocatable :: delta_OOVV(:,:,:,:)
double precision,allocatable :: OOOO(:,:,:,:)
double precision,allocatable :: OOOV(:,:,:,:)
double precision,allocatable :: OOVV(:,:,:,:)
double precision,allocatable :: OVOV(:,:,:,:)
double precision,allocatable :: OVVV(:,:,:,:)
double precision,allocatable :: VVVV(:,:,:,:)
double precision,allocatable :: X1(:,:,:,:)
double precision,allocatable :: X2(:,:)
double precision,allocatable :: X3(:,:)
double precision,allocatable :: X4(:,:,:,:)
double precision,allocatable :: u(:,:,:,:)
double precision,allocatable :: v(:,:,:,:)
double precision,allocatable :: r2(:,:,:,:)
double precision,allocatable :: t2(:,:,:,:)
integer :: n_diis,i,j,a,b
double precision :: rcond
double precision,allocatable :: error_diis(:,:)
double precision,allocatable :: t_diis(:,:)
! Hello world
write(*,*)
write(*,*)'**************************************'
write(*,*)'| BCCD calculation |'
write(*,*)'**************************************'
write(*,*)
! Spatial to spin orbitals
nBas = 2*nBasin
nC = 2*nCin
nO = 2*nOin
nV = 2*nVin
nR = 2*nRin
allocate(seHF(nBas),sERI(nBas,nBas,nBas,nBas))
call spatial_to_spin_MO_energy(nBasin,eHF,nBas,seHF)
call spatial_to_spin_ERI(nBasin,ERI,nBas,sERI)
! Antysymmetrize ERIs
allocate(dbERI(nBas,nBas,nBas,nBas))
call antisymmetrize_ERI(2,nBas,sERI,dbERI)
deallocate(sERI)
! Form energy denominator
allocate(eO(nO-nC),eV(nV-nR))
allocate(delta_OOVV(nO-nC,nO-nC,nV-nR,nV-nR))
eO(:) = seHF(nC+1:nO)
eV(:) = seHF(nO+1:nBas-nR)
call form_delta_OOVV(nC,nO,nV,nR,eO,eV,delta_OOVV)
deallocate(seHF)
! Create integral batches
allocate(OOOO(nO-nC,nO-nC,nO-nC,nO-nC),OOOV(nO-nC,nO-nC,nO-nC,nV-nR), &
OOVV(nO-nC,nO-nC,nV-nR,nV-nR),OVOV(nO-nC,nV-nR,nO-nC,nV-nR), &
OVVV(nO-nC,nV-nR,nV-nR,nV-nR),VVVV(nV-nR,nV-nR,nV-nR,nV-nR))
OOOO(:,:,:,:) = dbERI(nC+1:nO ,nC+1:nO ,nC+1:nO ,nC+1:nO )
OOOV(:,:,:,:) = dbERI(nC+1:nO ,nC+1:nO ,nC+1:nO ,nO+1:nBas-nR)
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)
OVVV(:,:,:,:) = dbERI(nC+1:nO ,nO+1:nBas-nR,nO+1:nBas-nR,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-nC,nO-nC,nV-nR,nV-nR))
t2(:,:,:,:) = -OOVV(:,:,:,:)/delta_OOVV(:,:,:,:)
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t2,EcMP2)
EcMP4 = 0d0
! Memory allocation for 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-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
n_diis = 0
t_diis(:,:) = 0d0
error_diis(:,:) = 0d0
!------------------------------------------------------------------------
! Main SCF loop
!------------------------------------------------------------------------
write(*,*)
write(*,*)'----------------------------------------------------'
write(*,*)'| BCCD calculation |'
write(*,*)'----------------------------------------------------'
write(*,'(1X,A1,1X,A3,1X,A1,1X,A16,1X,A1,1X,A10,1X,A1,1X,A10,1X,A1,1X)') &
'|','#','|','E(BCCD)','|','Ec(BCCD)','|','Conv','|'
write(*,*)'----------------------------------------------------'
do while(Conv > thresh .and. nSCF < maxSCF)
! Increment
nSCF = nSCF + 1
! Form linear array
call form_u(nC,nO,nV,nR,OOOO,VVVV,OVOV,t2,u)
! Form interemediate arrays
call form_X(nC,nO,nV,nR,OOVV,t2,X1,X2,X3,X4)
! Form quadratic array
call form_v(nC,nO,nV,nR,X1,X2,X3,X4,t2,v)
! Compute residual
r2(:,:,:,:) = OOVV(:,:,:,:) + delta_OOVV(:,:,:,:)*t2(:,:,:,:) + u(:,:,:,:) + v(:,:,:,:)
! Check convergence
Conv = maxval(abs(r2(:,:,:,:)))
! Update amplitudes
t2(:,:,:,:) = t2(:,:,:,:) - r2(:,:,:,:)/delta_OOVV(:,:,:,:)
! Compute correlation energy
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t2,EcCCD)
if(nSCF == 1) call MP3_correlation_energy(nC,nO,nV,nR,OOVV,t2,v,delta_OOVV,EcMP3)
! Dump results
ECCD = ERHF + EcCCD
! 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)
! Reset DIIS if required
if(abs(rcond) < 1d-15) n_diis = 0
write(*,'(1X,A1,1X,I3,1X,A1,1X,F16.10,1X,A1,1X,F10.6,1X,A1,1X,F10.6,1X,A1,1X)') &
'|',nSCF,'|',ECCD+ENuc,'|',EcCCD,'|',Conv,'|'
enddo
write(*,*)'----------------------------------------------------'
!------------------------------------------------------------------------
! End of SCF loop
!------------------------------------------------------------------------
! Did it actually converge?
if(nSCF == maxSCF) then
write(*,*)
write(*,*)'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
write(*,*)' Convergence failed '
write(*,*)'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
write(*,*)
stop
endif
! Moller-Plesset energies
write(*,*)
write(*,'(1X,A15,1X,F10.6)') 'Ec(MP2) = ',EcMP2
write(*,'(1X,A15,1X,F10.6)') 'Ec(MP3) = ',EcMP3
write(*,'(1X,A15,1X,F10.6)') 'Ec(MP4-SDQ) = ',EcMP4
write(*,*)
end subroutine

2
src/CC/EE_EOM_CCD_1h1p.f90
View File

@ -143,7 +143,7 @@ subroutine EE_EOM_CCD_1h1p(nC,nO,nV,nR,eO,eV,OOVV,OVVO,t)
call quick_sort(Om,order,nS)
call set_order(Z,order,nS,nS)
call print_excitation('EE-EOM-CCD ',3,nS,Om)
call print_excitation_energies('EE-EOM-CCD',3,nS,Om)
end if

12
src/CI/CIS.f90
View File

@ -55,8 +55,8 @@ subroutine CIS(singlet,triplet,doCIS_D,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF)
endif
call diagonalize_matrix(nS,A,Om)
call print_excitation('CIS ',ispin,nS,Om)
call print_transition_vectors_ph(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,transpose(A),transpose(A))
call print_excitation_energies('CIS',ispin,nS,Om)
call phLR_transition_vectors(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,transpose(A),transpose(A))
if(dump_trans) then
print*,'Singlet CIS transition vectors'
@ -67,7 +67,7 @@ subroutine CIS(singlet,triplet,doCIS_D,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF)
! Compute CIS(D) correction
maxS = min(maxS,nS)
if(doCIS_D) call D_correction(ispin,nBas,nC,nO,nV,nR,nS,maxS,eHF,ERI,Om(1:maxS),A(:,1:maxS))
if(doCIS_D) call CIS_D(ispin,nBas,nC,nO,nV,nR,nS,maxS,eHF,ERI,Om(1:maxS),A(:,1:maxS))
endif
@ -83,8 +83,8 @@ subroutine CIS(singlet,triplet,doCIS_D,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF)
endif
call diagonalize_matrix(nS,A,Om)
call print_excitation('CIS ',ispin,nS,Om)
call print_transition_vectors_ph(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,transpose(A),transpose(A))
call print_excitation_energies('CIS',ispin,nS,Om)
call phLR_transition_vectors(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,transpose(A),transpose(A))
if(dump_trans) then
print*,'Triplet CIS transition vectors'
@ -95,7 +95,7 @@ subroutine CIS(singlet,triplet,doCIS_D,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF)
! Compute CIS(D) correction
maxS = min(maxS,nS)
if(doCIS_D) call D_correction(ispin,nBas,nC,nO,nV,nR,nS,maxS,eHF,ERI,Om(1:maxS),A(:,1:maxS))
if(doCIS_D) call CIS_D(ispin,nBas,nC,nO,nV,nR,nS,maxS,eHF,ERI,Om(1:maxS),A(:,1:maxS))
endif

3
src/CI/D_correction.f90 → src/CI/CIS_D.f90
View File

@ -1,4 +1,4 @@
subroutine D_correction(ispin,nBasin,nCin,nOin,nVin,nRin,nSin,maxS,eHF,ERI,w,X)
subroutine CIS_D(ispin,nBasin,nCin,nOin,nVin,nRin,nSin,maxS,eHF,ERI,w,X)
! Compute the D correction of CIS(D)
@ -274,4 +274,5 @@ subroutine D_correction(ispin,nBasin,nCin,nOin,nVin,nRin,nSin,maxS,eHF,ERI,w,X)
!------------------------------------------------------------------------
! End of loop over single excitations
!------------------------------------------------------------------------
end subroutine

28
src/CI/UCIS.f90
View File

@ -34,11 +34,11 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E
integer :: nS_aa,nS_bb,nS_sc
double precision,allocatable :: A_sc(:,:)
double precision,allocatable :: Omega_sc(:)
double precision,allocatable :: Om_sc(:)
integer :: nS_ab,nS_ba,nS_sf
double precision,allocatable :: A_sf(:,:)
double precision,allocatable :: Omega_sf(:)
double precision,allocatable :: Om_sf(:)
! Hello world
@ -66,7 +66,7 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E
nS_bb = nS(2)
nS_sc = nS_aa + nS_bb
allocate(A_sc(nS_sc,nS_sc),Omega_sc(nS_sc))
allocate(A_sc(nS_sc,nS_sc),Om_sc(nS_sc))
call phULR_A(ispin,.false.,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,eHF,ERI_aaaa,ERI_aabb,ERI_bbbb,A_sc)
@ -76,11 +76,11 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E
write(*,*)
endif
call diagonalize_matrix(nS_sc,A_sc,Omega_sc)
call diagonalize_matrix(nS_sc,A_sc,Om_sc)
A_sc(:,:) = transpose(A_sc)
call print_excitation('UCIS ',5,nS_sc,Omega_sc)
call print_unrestricted_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_aa,nS_bb,nS_sc,dipole_int_aa,dipole_int_bb, &
cHF,S,Omega_sc,A_sc,A_sc)
call print_excitation_energies('UCIS',5,nS_sc,Om_sc)
call phULR_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_aa,nS_bb,nS_sc,dipole_int_aa,dipole_int_bb, &
cHF,S,Om_sc,A_sc,A_sc)
if(dump_trans) then
print*,'Spin-conserved CIS transition vectors'
@ -88,7 +88,7 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E
write(*,*)
endif
deallocate(A_sc,Omega_sc)
deallocate(A_sc,Om_sc)
endif
@ -106,7 +106,7 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E
nS_ba = (nO(2) - nC(2))*(nV(1) - nR(1))
nS_sf = nS_ab + nS_ba
allocate(A_sf(nS_sf,nS_sf),Omega_sf(nS_sf))
allocate(A_sf(nS_sf,nS_sf),Om_sf(nS_sf))
call phULR_A(ispin,.false.,nBas,nC,nO,nV,nR,nS_ab,nS_ba,nS_sf,lambda,eHF,ERI_aaaa,ERI_aabb,ERI_bbbb,A_sf)
@ -116,11 +116,11 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E
write(*,*)
endif
call diagonalize_matrix(nS_sf,A_sf,Omega_sf)
call diagonalize_matrix(nS_sf,A_sf,Om_sf)
A_sf(:,:) = transpose(A_sf)
call print_excitation('UCIS ',6,nS_sf,Omega_sf)
call print_unrestricted_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_ab,nS_ba,nS_sf,dipole_int_aa,dipole_int_bb, &
cHF,S,Omega_sf,A_sf,A_sf)
call print_excitation_energies('UCIS ',6,nS_sf,Om_sf)
call phULR_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_ab,nS_ba,nS_sf,dipole_int_aa,dipole_int_bb, &
cHF,S,Om_sf,A_sf,A_sf)
if(dump_trans) then
print*,'Spin-flip CIS transition vectors'
@ -128,7 +128,7 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E
write(*,*)
endif
deallocate(A_sf,Omega_sf)
deallocate(A_sf,Om_sf)
endif

8
src/GF/GF2_phBSE2.f90
View File

@ -68,8 +68,8 @@ subroutine GF2_phBSE2(TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,
! Compute phBSE2@GF2 excitation energies
call phLR(TDA,nS,A_sta,B_sta,EcBSE(ispin),OmBSE,XpY,XmY)
call print_excitation('phBSE2@GF2 ',ispin,nS,OmBSE)
call print_transition_vectors_ph(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY,XmY)
call print_excitation_energies('phBSE2@GF2',ispin,nS,OmBSE)
call phLR_transition_vectors(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY,XmY)
! Compute dynamic correction for BSE via perturbation theory
@ -101,8 +101,8 @@ subroutine GF2_phBSE2(TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,
! Compute phBSE2@GF2 excitation energies
call phLR(TDA,nS,A_sta,B_sta,EcBSE(ispin),OmBSE,XpY,XmY)
call print_excitation('phBSE2@GF2 ',ispin,nS,OmBSE)
call print_transition_vectors_ph(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY,XmY)
call print_excitation_energies('phBSE2@GF2',ispin,nS,OmBSE)
call phLR_transition_vectors(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY,XmY)
! Compute dynamic correction for BSE via perturbation theory

2
src/GF/GF2_phBSE2_dynamic_perturbation.f90
View File

@ -57,7 +57,7 @@ subroutine GF2_phBSE2_dynamic_perturbation(dTDA,ispin,eta,nBas,nC,nO,nV,nR,nS,ER
! Print main components of transition vectors
call print_transition_vectors_ph(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY,XmY)
call phLR_transition_vectors(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY,XmY)
gapGF = eGF(nO+1) - eGF(nO)

4
src/GF/GF2_ppBSE2.f90
View File

@ -93,7 +93,7 @@ subroutine GF2_ppBSE2(TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,ERI,dip
call ppLR(TDA,nOO,nVV,Bpp,Cpp,Dpp,Om1,X1,Y1,Om2,X2,Y2,EcBSE(ispin))
call print_transition_vectors_pp(.true.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
call ppLR_transition_vectors(.true.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
!----------------------------------------------------!
! Compute the dynamical screening at the ppBSE level !
@ -146,7 +146,7 @@ subroutine GF2_ppBSE2(TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,ERI,dip
call ppLR(TDA,nOO,nVV,Bpp,Cpp,Dpp,Om1,X1,Y1,Om2,X2,Y2,EcBSE(ispin))
call print_transition_vectors_pp(.false.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
call ppLR_transition_vectors(.false.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
!----------------------------------------------------!
! Compute the dynamical screening at the ppBSE level !

2
src/GT/G0T0eh.f90
View File

@ -114,7 +114,7 @@ subroutine G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE,
call phLR(TDA_T,nS,Aph,Bph,EcRPA,Om,XpY,XmY)
if(print_T) call print_excitation('RPA@HF ',ispin,nS,Om)
if(print_T) call print_excitation_energies('RPA@HF ',ispin,nS,Om)
!--------------------------!
! Compute spectral weights !

8
src/GT/G0T0pp.f90
View File

@ -109,8 +109,8 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,dopp
deallocate(Bpp,Cpp,Dpp)
call print_excitation('pp-RPA (N+2)',iblock,nVVs,Om1s(:))
call print_excitation('pp-RPA (N-2)',iblock,nOOs,Om2s(:))
call print_excitation_energies('pp-RPA (N+2)',iblock,nVVs,Om1s(:))
call print_excitation_energies('pp-RPA (N-2)',iblock,nOOs,Om2s(:))
!----------------------------------------------
! alpha-alpha block
@ -131,8 +131,8 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,dopp
deallocate(Bpp,Cpp,Dpp)
call print_excitation('ppRPA (N+2) ',iblock,nVVt,Om1t)
call print_excitation('ppRPA (N-2) ',iblock,nOOt,Om2t)
call print_excitation_energies('ppRPA (N+2)',iblock,nVVt,Om1t)
call print_excitation_energies('ppRPA (N-2)',iblock,nOOt,Om2t)
!----------------------------------------------
! Compute T-matrix version of the self-energy

8
src/GT/GTpp_phBSE.f90
View File

@ -139,8 +139,8 @@ subroutine GTpp_phBSE(TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,n
call phLR(TDA,nS,Aph,Bph,EcBSE(ispin),OmBSE,XpY_BSE,XmY_BSE)
call print_excitation('phBSE@GTpp ',ispin,nS,OmBSE)
call print_transition_vectors_ph(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY_BSE,XmY_BSE)
call print_excitation_energies('phBSE@GTpp',ispin,nS,OmBSE)
call phLR_transition_vectors(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY_BSE,XmY_BSE)
! Compute dynamic correction for BSE via renormalized perturbation theory
@ -169,8 +169,8 @@ subroutine GTpp_phBSE(TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,n
call phLR(TDA,nS,Aph,Bph,EcBSE(ispin),OmBSE,XpY_BSE,XmY_BSE)
call print_excitation('phBSE@GTpp ',ispin,nS,OmBSE)
call print_transition_vectors_ph(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY_BSE,XmY_BSE)
call print_excitation_energies('phBSE@GTpp',ispin,nS,OmBSE)
call phLR_transition_vectors(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY_BSE,XmY_BSE)
! Compute dynamic correction for BSE via renormalized perturbation theory

4
src/GT/GTpp_ppBSE.f90
View File

@ -155,7 +155,7 @@ subroutine GTpp_ppBSE(TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,n
call ppLR(TDA,nOOs,nVVs,Bpp,Cpp,Dpp,Om1s,X1s,Y1s,Om2s,X2s,Y2s,EcBSE(ispin))
call print_transition_vectors_pp(.true.,nBas,nC,nO,nV,nR,nOOs,nVVs,dipole_int,Om1s,X1s,Y1s,Om2s,X2s,Y2s)
call ppLR_transition_vectors(.true.,nBas,nC,nO,nV,nR,nOOs,nVVs,dipole_int,Om1s,X1s,Y1s,Om2s,X2s,Y2s)
deallocate(Om1s,X1s,Y1s,Om2s,X2s,Y2s,TBab,TCab,TDab,TBaa,TCaa,TDaa,Bpp,Cpp,Dpp)
@ -240,7 +240,7 @@ subroutine GTpp_ppBSE(TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,n
call ppLR(TDA,nOOs,nVVs,Bpp,Cpp,Dpp,Om1s,X1s,Y1s,Om2s,X2s,Y2s,EcBSE(ispin))
call print_transition_vectors_pp(.false.,nBas,nC,nO,nV,nR,nOOt,nVVt,dipole_int,Om1t,X1t,Y1t,Om2t,X2t,Y2t)
call ppLR_transition_vectors(.false.,nBas,nC,nO,nV,nR,nOOt,nVVt,dipole_int,Om1t,X1t,Y1t,Om2t,X2t,Y2t)
deallocate(Om1t,X1t,Y1t,Om2t,X2t,Y2t,TBab,TCab,TDab,TBaa,TCaa,TDaa,Bpp,Cpp,Dpp)

12
src/GT/UG0T0pp.f90
View File

@ -120,8 +120,8 @@ subroutine UG0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA, &
! EcRPA(ispin) = 1d0*EcRPA(ispin)
call print_excitation('pp-RPA (N+2)',iblock,nPab,Om1ab(:))
call print_excitation('pp-RPA (N-2)',iblock,nHab,Om2ab(:))
call print_excitation_energies('ppRPA@HF (N+2)',iblock,nPab,Om1ab(:))
call print_excitation_energies('ppRPA@HF (N-2)',iblock,nHab,Om2ab(:))
!----------------------------------------------
! alpha-alpha block
@ -138,8 +138,8 @@ subroutine UG0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA, &
! EcRPA(ispin) = 2d0*EcRPA(ispin)
! EcRPA(ispin) = 3d0*EcRPA(ispin)
call print_excitation('pp-RPA (N+2)',iblock,nPaa,Om1aa(:))
call print_excitation('pp-RPA (N-2)',iblock,nHaa,Om2aa(:))
call print_excitation_energies('ppRPA@HF (N+2)',iblock,nPaa,Om1aa(:))
call print_excitation_energies('ppRPA@HF (N-2)',iblock,nHaa,Om2aa(:))
!----------------------------------------------
! beta-beta block
@ -156,8 +156,8 @@ subroutine UG0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA, &
! EcRPA(ispin) = 2d0*EcRPA(ispin)
! EcRPA(ispin) = 3d0*EcRPA(ispin)
call print_excitation('pp-RPA (N+2)',iblock,nPbb,Om1bb(:))
call print_excitation('pp-RPA (N-2)',iblock,nHbb,Om2bb(:))
call print_excitation_energies('ppRPA@HF (N+2)',iblock,nPbb,Om1bb(:))
call print_excitation_energies('ppRPA@HF (N-2)',iblock,nHbb,Om2bb(:))
!----------------------------------------------
! Compute T-matrix version of the self-energy

2
src/GT/qsGTeh.f90
View File

@ -177,7 +177,7 @@ subroutine qsGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE2,
call phLR(TDA_T,nS,Aph,Bph,EcRPA,Om,XpY,XmY)
if(print_T) call print_excitation('RPA@qsGTeh ',ispin,nS,Om)
if(print_T) call print_excitation_energies('phRPA@qsGTeh',ispin,nS,Om)
! Compute correlation part of the self-energy

2
src/GW/G0W0.f90
View File

@ -104,7 +104,7 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dT
call phLR(TDA_W,nS,Aph,Bph,EcRPA,Om,XpY,XmY)
if(print_W) call print_excitation('RPA@HF ',ispin,nS,Om)
if(print_W) call print_excitation_energies('RPA@HF ',ispin,nS,Om)
!--------------------------!
! Compute spectral weights !

8
src/GW/GW_phBSE.f90
View File

@ -117,8 +117,8 @@ subroutine GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,
call phLR(TDA,nS,Aph,Bph,EcBSE(ispin),OmBSE,XpY_BSE,XmY_BSE)
call print_excitation('phBSE@GW ',ispin,nS,OmBSE)
call print_transition_vectors_ph(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY_BSE,XmY_BSE)
call print_excitation_energies('phBSE@GW',ispin,nS,OmBSE)
call phLR_transition_vectors(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY_BSE,XmY_BSE)
!----------------------------------------------------!
! Compute the dynamical screening at the phBSE level !
@ -149,8 +149,8 @@ subroutine GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,
call phLR(TDA,nS,Aph,Bph,EcBSE(ispin),OmBSE,XpY_BSE,XmY_BSE)
call print_excitation('phBSE@GW ',ispin,nS,OmBSE)
call print_transition_vectors_ph(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY_BSE,XmY_BSE)
call print_excitation_energies('phBSE@GW',ispin,nS,OmBSE)
call phLR_transition_vectors(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY_BSE,XmY_BSE)
!-------------------------------------------------
! Compute the dynamical screening at the BSE level

4
src/GW/GW_ppBSE.f90
View File

@ -122,7 +122,7 @@ subroutine GW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,
call ppLR(TDA,nOO,nVV,Bpp,Cpp,Dpp,Om1,X1,Y1,Om2,X2,Y2,EcBSE(ispin))
call print_transition_vectors_pp(.true.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
call ppLR_transition_vectors(.true.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
!----------------------------------------------------!
! Compute the dynamical screening at the ppBSE level !
@ -175,7 +175,7 @@ subroutine GW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,
call ppLR(TDA,nOO,nVV,Bpp,Cpp,Dpp,Om1,X1,Y1,Om2,X2,Y2,EcBSE(ispin))
call print_transition_vectors_pp(.false.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
call ppLR_transition_vectors(.false.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
!----------------------------------------------------!
! Compute the dynamical screening at the ppBSE level !

4
src/GW/SRG_qsGW.f90
View File

@ -141,8 +141,6 @@ subroutine SRG_qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE
error_diis(:,:) = 0d0
rcond = 0d0
print*,max_diis
!------------------------------------------------------------------------
! Main loop
!------------------------------------------------------------------------
@ -183,7 +181,7 @@ subroutine SRG_qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE
tlr = tlr + tlr2 -tlr1
if(print_W) call print_excitation('RPA@qsGW ',ispin,nS,OmRPA)
if(print_W) call print_excitation_energies('phRPA@SRG-qsGW',ispin,nS,OmRPA)
! Compute correlation part of the self-energy

2
src/GW/UG0W0.f90
View File

@ -112,7 +112,7 @@ subroutine UG0W0(doACFDT,exchange_kernel,doXBS,BSE,TDA_W,TDA,dBSE,dTDA,spin_cons
call phULR(ispin,.true.,TDA_W,.false.,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0, &
eHF,ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
if(print_W) call print_excitation('RPA@UHF ',5,nS_sc,OmRPA)
if(print_W) call print_excitation_energies('phRPA@UHF',5,nS_sc,OmRPA)
!----------------------!
! Excitation densities !

12
src/GW/UGW_phBSE.f90
View File

@ -100,9 +100,9 @@ subroutine UGW_phBSE(TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta, &
call phULR(ispin,.true.,TDA,.true.,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0, &
eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,EcBSE(ispin),OmBSE_sc,XpY_BSE_sc,XmY_BSE_sc)
call print_excitation('BSE@UGW ',5,nS_sc,OmBSE_sc)
call print_unrestricted_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_aa,nS_bb,nS_sc,dipole_int_aa,dipole_int_bb, &
cW,S,OmBSE_sc,XpY_BSE_sc,XmY_BSE_sc)
call print_excitation_energies('phBSE@UGW',5,nS_sc,OmBSE_sc)
call phULR_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_aa,nS_bb,nS_sc,dipole_int_aa,dipole_int_bb, &
cW,S,OmBSE_sc,XpY_BSE_sc,XmY_BSE_sc)
!-------------------------------------------------
! Compute the dynamical screening at the BSE level
@ -136,9 +136,9 @@ subroutine UGW_phBSE(TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta, &
eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,EcBSE(ispin), &
OmBSE_sf,XpY_BSE_sf,XmY_BSE_sf)
call print_excitation('BSE@UGW ',6,nS_sf,OmBSE_sf)
call print_unrestricted_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_ab,nS_ba,nS_sf,dipole_int_aa,dipole_int_bb, &
cW,S,OmBSE_sf,XpY_BSE_sf,XmY_BSE_sf)
call print_excitation_energies('phBSE@UGW',6,nS_sf,OmBSE_sf)
call phULR_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_ab,nS_ba,nS_sf,dipole_int_aa,dipole_int_bb, &
cW,S,OmBSE_sf,XpY_BSE_sf,XmY_BSE_sf)
!-------------------------------------------------
! Compute the dynamical screening at the BSE level

2
src/GW/qsGW.f90
View File

@ -176,7 +176,7 @@ subroutine qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dop
if(.not.TDA_W) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI_MO,Bph)
call phLR(TDA_W,nS,Aph,Bph,EcRPA,Om,XpY,XmY)
if(print_W) call print_excitation('phRPA@qsGW ',ispin,nS,Om)
if(print_W) call print_excitation_energies('phRPA@qsGW',ispin,nS,Om)
! Compute correlation part of the self-energy

2
src/LR/oscillator_strength_ph.f90 → src/LR/phLR_oscillator_strength.f90
View File

@ -1,4 +1,4 @@
subroutine oscillator_strength_ph(nBas,nC,nO,nV,nR,nS,maxS,dipole_int,Om,XpY,XmY,os)
subroutine phLR_oscillator_strength(nBas,nC,nO,nV,nR,nS,maxS,dipole_int,Om,XpY,XmY,os)
! Compute linear response

4
src/LR/print_transition_vectors_ph.f90 → src/LR/phLR_transition_vectors.f90
View File

@ -1,4 +1,4 @@
subroutine print_transition_vectors_ph(spin_allowed,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
subroutine phLR_transition_vectors(spin_allowed,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
! Print transition vectors for linear response calculation
@ -37,7 +37,7 @@ subroutine print_transition_vectors_ph(spin_allowed,nBas,nC,nO,nV,nR,nS,dipole_i
! Compute oscillator strengths
os(:) = 0d0
if(spin_allowed) call oscillator_strength_ph(nBas,nC,nO,nV,nR,nS,maxS,dipole_int,Om,XpY,XmY,os)
if(spin_allowed) call phLR_oscillator_strength(nBas,nC,nO,nV,nR,nS,maxS,dipole_int,Om,XpY,XmY,os)
! Print details about excitations

2
src/LR/unrestricted_oscillator_strength.f90 → src/LR/phULR_oscillator_strength.f90
View File

@ -1,4 +1,4 @@
subroutine unrestricted_oscillator_strength(nBas,nC,nO,nV,nR,nS,nSa,nSb,nSt,maxS,dipole_int_aa,dipole_int_bb,Omega,XpY,XmY,os)
subroutine phULR_oscillator_strength(nBas,nC,nO,nV,nR,nS,nSa,nSb,nSt,maxS,dipole_int_aa,dipole_int_bb,Omega,XpY,XmY,os)
! Compute linear response

7
src/LR/print_unrestricted_transition_vectors.f90 → src/LR/phULR_transition_vectors.f90
View File

@ -1,5 +1,4 @@
subroutine print_unrestricted_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nSa,nSb,nSt,dipole_int_aa,dipole_int_bb, &
c,S,Omega,XpY,XmY)
subroutine phULR_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nSa,nSb,nSt,dipole_int_aa,dipole_int_bb,c,S,Omega,XpY,XmY)
! Print transition vectors for linear response calculation
@ -44,8 +43,8 @@ subroutine print_unrestricted_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nSa,n
! Compute oscillator strengths
os(:) = 0d0
if(ispin == 1) call unrestricted_oscillator_strength(nBas,nC,nO,nV,nR,nS,nSa,nSb,nSt,maxS, &
dipole_int_aa,dipole_int_bb,Omega,XpY,XmY,os)
if(ispin == 1) call phULR_oscillator_strength(nBas,nC,nO,nV,nR,nS,nSa,nSb,nSt,maxS, &
dipole_int_aa,dipole_int_bb,Omega,XpY,XmY,os)
! Compute <S**2>

2
src/LR/oscillator_strength_pp.f90 → src/LR/ppLR_oscillator_strength.f90
View File

@ -1,4 +1,4 @@
subroutine oscillator_strength_pp(nBas,nC,nO,nV,nR,nOO,nVV,maxOO,maxVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2,os1,os2)
subroutine ppLR_oscillator_strength(nBas,nC,nO,nV,nR,nOO,nVV,maxOO,maxVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2,os1,os2)
! Compute linear response

4
src/LR/print_transition_vectors_pp.f90 → src/LR/ppLR_transition_vectors.f90
View File

@ -1,4 +1,4 @@
subroutine print_transition_vectors_pp(spin_allowed,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
subroutine ppLR_transition_vectors(spin_allowed,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
! Print transition vectors for p-p linear response calculation
@ -47,7 +47,7 @@ subroutine print_transition_vectors_pp(spin_allowed,nBas,nC,nO,nV,nR,nOO,nVV,dip
os2(:) = 0d0
if(spin_allowed) &
call oscillator_strength_pp(nBas,nC,nO,nV,nR,nOO,nVV,maxOO,maxVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2,os1,os2)
call ppLR_oscillator_strength(nBas,nC,nO,nV,nR,nOO,nVV,maxOO,maxVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2,os1,os2)
!-----------------------------------------------!
! Print details about excitations for pp sector !

25
src/LR/print_excitation.f90 → src/LR/print_excitation_energies.f90
View File

@ -1,4 +1,4 @@
subroutine print_excitation(method,ispin,nS,Omega)
subroutine print_excitation_energies(method,ispin,nS,Om)
! Print excitation energies for a given spin manifold
@ -7,14 +7,15 @@ subroutine print_excitation(method,ispin,nS,Omega)
! Input variables
character*12,intent(in) :: method
integer,intent(in) :: ispin,nS
double precision,intent(in) :: Omega(nS)
character(len=20),intent(in) :: method
integer,intent(in) :: ispin
integer,intent(in) :: nS
double precision,intent(in) :: Om(nS)
! Local variables
character*14 :: spin_manifold
integer,parameter :: maxS = 50
character(len=20) :: spin_manifold
integer,parameter :: maxS = 20
integer :: ia
if(ispin == 1) spin_manifold = 'singlet'
@ -26,19 +27,19 @@ subroutine print_excitation(method,ispin,nS,Omega)
if(ispin == 7) spin_manifold = 'beta-beta'
write(*,*)
write(*,*)'-------------------------------------------------------------'
write(*,'(1X,A14,A14,A14,A9)') method,' calculation: ',spin_manifold,' manifold'
write(*,*)'-------------------------------------------------------------'
write(*,*)'-------------------------------------------------------------------------------'
write(*,'(1X,A20,A20,A20,A9)') trim(method),' calculation: ',trim(spin_manifold),' manifold'
write(*,*)'-------------------------------------------------------------------------------'
write(*,'(1X,A1,1X,A5,1X,A1,1X,A23,1X,A1,1X,A23,1X,A1,1X)') &
'|','State','|',' Excitation energy (au) ','|',' Excitation energy (eV) ','|'
write(*,*)'-------------------------------------------------------------'
write(*,*)'-------------------------------------------------------------------------------'
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,'|'
'|',ia,'|',Om(ia),'|',Om(ia)*HaToeV,'|'
enddo
write(*,*)'-------------------------------------------------------------'
write(*,*)'-------------------------------------------------------------------------------'
write(*,*)
end subroutine

8
src/RPA/crRPA.f90
View File

@ -75,8 +75,8 @@ subroutine crRPA(TDA,doACFDT,exchange_kernel,singlet,triplet,nBas,nC,nO,nV,nR,nS
if(.not.TDA) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,-1d0,ERI,Bph)
call phLR(TDA,nS,Aph,Bph,EcTr(ispin),Om,XpY,XmY)
call print_excitation('crRPA@HF ',ispin,nS,Om)
call print_transition_vectors_ph(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
call print_excitation_energies('crRPA@HF',ispin,nS,Om)
call phLR_transition_vectors(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
endif
@ -90,8 +90,8 @@ subroutine crRPA(TDA,doACFDT,exchange_kernel,singlet,triplet,nBas,nC,nO,nV,nR,nS
if(.not.TDA) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,-1d0,ERI,Bph)
call phLR(TDA,nS,Aph,Bph,EcTr(ispin),Om,XpY,XmY)
call print_excitation('crRPA@HF ',ispin,nS,Om)
call print_transition_vectors_ph(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
call print_excitation_energies('crRPA@HF',ispin,nS,Om)
call phLR_transition_vectors(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
endif

8
src/RPA/phRPA.f90
View File

@ -75,8 +75,8 @@ subroutine phRPA(TDA,doACFDT,exchange_kernel,singlet,triplet,nBas,nC,nO,nV,nR,nS
if(.not.TDA) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph)
call phLR(TDA,nS,Aph,Bph,EcTr(ispin),Om,XpY,XmY)
call print_excitation('phRPA@HF ',ispin,nS,Om)
call print_transition_vectors_ph(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
call print_excitation_energies('phRPA@HF',ispin,nS,Om)
call phLR_transition_vectors(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
endif
@ -90,8 +90,8 @@ subroutine phRPA(TDA,doACFDT,exchange_kernel,singlet,triplet,nBas,nC,nO,nV,nR,nS
if(.not.TDA) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph)
call phLR(TDA,nS,Aph,Bph,EcTr(ispin),Om,XpY,XmY)
call print_excitation('phRPA@HF ',ispin,nS,Om)
call print_transition_vectors_ph(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
call print_excitation_energies('phRPA@HF ',ispin,nS,Om)
call phLR_transition_vectors(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
endif

8
src/RPA/phRPAx.f90
View File

@ -76,8 +76,8 @@ subroutine phRPAx(TDA,doACFDT,exchange_kernel,singlet,triplet,nBas,nC,nO,nV,nR,n
if(.not.TDA) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph)
call phLR(TDA,nS,Aph,Bph,EcTr(ispin),Om,XpY,XmY)
call print_excitation('phRPAx@HF ',ispin,nS,Om)
call print_transition_vectors_ph(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
call print_excitation_energies('phRPAx@HF ',ispin,nS,Om)
call phLR_transition_vectors(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
endif
@ -91,8 +91,8 @@ subroutine phRPAx(TDA,doACFDT,exchange_kernel,singlet,triplet,nBas,nC,nO,nV,nR,n
if(.not.TDA) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph)
call phLR(TDA,nS,Aph,Bph,EcTr(ispin),Om,XpY,XmY)
call print_excitation('phRPAx@HF ',ispin,nS,Om)
call print_transition_vectors_ph(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
call print_excitation_energies('phRPAx@HF ',ispin,nS,Om)
call phLR_transition_vectors(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,Om,XpY,XmY)
endif

28
src/RPA/phURPA.f90
View File

@ -36,12 +36,12 @@ subroutine phURPA(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,nBas,nC,n
integer :: ispin
integer :: nS_aa,nS_bb,nS_sc
double precision,allocatable :: Omega_sc(:)
double precision,allocatable :: Om_sc(:)
double precision,allocatable :: XpY_sc(:,:)
double precision,allocatable :: XmY_sc(:,:)
integer :: nS_ab,nS_ba,nS_sf
double precision,allocatable :: Omega_sf(:)
double precision,allocatable :: Om_sf(:)
double precision,allocatable :: XpY_sf(:,:)
double precision,allocatable :: XmY_sf(:,:)
@ -80,15 +80,15 @@ subroutine phURPA(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,nBas,nC,n
nS_bb = nS(2)
nS_sc = nS_aa + nS_bb
allocate(Omega_sc(nS_sc),XpY_sc(nS_sc,nS_sc),XmY_sc(nS_sc,nS_sc))
allocate(Om_sc(nS_sc),XpY_sc(nS_sc,nS_sc),XmY_sc(nS_sc,nS_sc))
call phULR(ispin,.true.,TDA,.false.,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0,e, &
ERI_aaaa,ERI_aabb,ERI_bbbb,Omega_sc,rho_sc,EcRPA(ispin),Omega_sc,XpY_sc,XmY_sc)
call print_excitation('URPA ',5,nS_sc,Omega_sc)
call print_unrestricted_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_aa,nS_bb,nS_sc,dipole_int_aa,dipole_int_bb, &
c,S,Omega_sc,XpY_sc,XmY_sc)
ERI_aaaa,ERI_aabb,ERI_bbbb,Om_sc,rho_sc,EcRPA(ispin),Om_sc,XpY_sc,XmY_sc)
call print_excitation_energies('phURPA@HF',5,nS_sc,Om_sc)
call phULR_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_aa,nS_bb,nS_sc,dipole_int_aa,dipole_int_bb, &
c,S,Om_sc,XpY_sc,XmY_sc)
deallocate(Omega_sc,XpY_sc,XmY_sc)
deallocate(Om_sc,XpY_sc,XmY_sc)
endif
@ -104,15 +104,15 @@ subroutine phURPA(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,nBas,nC,n
nS_ba = (nO(2) - nC(2))*(nV(1) - nR(1))
nS_sf = nS_ab + nS_ba
allocate(Omega_sf(nS_sf),XpY_sf(nS_sf,nS_sf),XmY_sf(nS_sf,nS_sf))
allocate(Om_sf(nS_sf),XpY_sf(nS_sf,nS_sf),XmY_sf(nS_sf,nS_sf))
call phULR(ispin,.true.,TDA,.false.,nBas,nC,nO,nV,nR,nS_ab,nS_ba,nS_sf,nS_sf,1d0,e, &
ERI_aaaa,ERI_aabb,ERI_bbbb,Omega_sf,rho_sf,EcRPA(ispin),Omega_sf,XpY_sf,XmY_sf)
call print_excitation('URPA ',6,nS_sf,Omega_sf)
call print_unrestricted_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_ab,nS_ba,nS_sf,dipole_int_aa,dipole_int_bb, &
c,S,Omega_sf,XpY_sf,XmY_sf)
ERI_aaaa,ERI_aabb,ERI_bbbb,Om_sf,rho_sf,EcRPA(ispin),Om_sf,XpY_sf,XmY_sf)
call print_excitation_energies('phURPA@HF',6,nS_sf,Om_sf)
call phULR_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_ab,nS_ba,nS_sf,dipole_int_aa,dipole_int_bb, &
c,S,Om_sf,XpY_sf,XmY_sf)
deallocate(Omega_sf,XpY_sf,XmY_sf)
deallocate(Om_sf,XpY_sf,XmY_sf)
endif

12
src/RPA/phURPAx.f90
View File

@ -85,9 +85,9 @@ subroutine phURPAx(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,nBas,nC,
call phULR(ispin,.false.,TDA,.false.,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0,e, &
ERI_aaaa,ERI_aabb,ERI_bbbb,Omega_sc,rho_sc,EcRPA(ispin),Omega_sc,XpY_sc,XmY_sc)
call print_excitation('URPAx ',5,nS_sc,Omega_sc)
call print_unrestricted_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_aa,nS_bb,nS_sc,dipole_int_aa,dipole_int_bb, &
c,S,Omega_sc,XpY_sc,XmY_sc)
call print_excitation_energies('phURPAx@HF',5,nS_sc,Omega_sc)
call phULR_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_aa,nS_bb,nS_sc,dipole_int_aa,dipole_int_bb, &
c,S,Omega_sc,XpY_sc,XmY_sc)
deallocate(Omega_sc,XpY_sc,XmY_sc)
@ -109,9 +109,9 @@ subroutine phURPAx(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,nBas,nC,
call phULR(ispin,.false.,TDA,.false.,nBas,nC,nO,nV,nR,nS_ab,nS_ba,nS_sf,nS_sf,1d0,e, &
ERI_aaaa,ERI_aabb,ERI_bbbb,Omega_sf,rho_sf,EcRPA(ispin),Omega_sf,XpY_sf,XmY_sf)
call print_excitation('URPAx ',6,nS_sf,Omega_sf)
call print_unrestricted_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_ab,nS_ba,nS_sf,dipole_int_aa,dipole_int_bb, &
c,S,Omega_sf,XpY_sf,XmY_sf)
call print_excitation_energies('phURPAx@HF',6,nS_sf,Omega_sf)
call phULR_transition_vectors(ispin,nBas,nC,nO,nV,nR,nS,nS_ab,nS_ba,nS_sf,dipole_int_aa,dipole_int_bb, &
c,S,Omega_sf,XpY_sf,XmY_sf)
deallocate(Omega_sf,XpY_sf,XmY_sf)

8
src/RPA/ppRPA.f90
View File

@ -78,8 +78,8 @@ subroutine ppRPA(TDA,doACFDT,singlet,triplet,nBas,nC,nO,nV,nR,ENuc,EHF,ERI,dipol
! call print_transition_vectors_pp(.true.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
call print_excitation('ppRPA (N+2) ',ispin,nVV,Om1)
call print_excitation('ppRPA (N-2) ',ispin,nOO,Om2)
call print_excitation_energies('ppRPA (N+2) ',ispin,nVV,Om1)
call print_excitation_energies('ppRPA (N-2) ',ispin,nOO,Om2)
deallocate(Om1,X1,Y1,Om2,X2,Y2,Bpp,Cpp,Dpp)
@ -110,8 +110,8 @@ subroutine ppRPA(TDA,doACFDT,singlet,triplet,nBas,nC,nO,nV,nR,ENuc,EHF,ERI,dipol
! call print_transition_vectors_pp(.false.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
call print_excitation('ppRPA (N+2) ',ispin,nVV,Om1)
call print_excitation('ppRPA (N-2) ',ispin,nOO,Om2)
call print_excitation_energies('ppRPA (N+2) ',ispin,nVV,Om1)
call print_excitation_energies('ppRPA (N-2) ',ispin,nOO,Om2)
deallocate(Om1,X1,Y1,Om2,X2,Y2,Bpp,Cpp,Dpp)

12
src/RPA/ppURPA.f90
View File

@ -72,8 +72,8 @@ subroutine ppURPA(TDA,doACFDT,spin_conserved,spin_flip,nBas,nC,nO,nV,nR,ENuc,EUH
ERI_aabb,ERI_bbbb,Om1sc,X1sc,Y1sc, &
Om2sc,X2sc,Y2sc,Ec_ppURPA(ispin))
call print_excitation('pp-RPA (N+2)',iblock,nP_sc,Om1sc)
call print_excitation('pp-RPA (N-2)',iblock,nH_sc,Om2sc)
call print_excitation_energies('ppRPA@HF (N+2)',iblock,nP_sc,Om1sc)
call print_excitation_energies('ppRPA@HF (N-2)',iblock,nH_sc,Om2sc)
!alpha-alpha block
@ -98,8 +98,8 @@ subroutine ppURPA(TDA,doACFDT,spin_conserved,spin_flip,nBas,nC,nO,nV,nR,ENuc,EUH
ERI_aabb,ERI_bbbb,Om1sf,X1sf,Y1sf, &
Om2sf,X2sf,Y2sf,Ec_ppURPA(ispin))
call print_excitation('pp-RPA (N+2)',iblock,nP_sf,Om1sf)
call print_excitation('pp-RPA (N-2)',iblock,nH_sf,Om2sf)
call print_excitation_energies('ppRPA@HF (N+2)',iblock,nP_sf,Om1sf)
call print_excitation_energies('ppRPA@HF (N-2)',iblock,nH_sf,Om2sf)
deallocate(Om1sf,X1sf,Y1sf,Om2sf,X2sf,Y2sf)
@ -118,8 +118,8 @@ subroutine ppURPA(TDA,doACFDT,spin_conserved,spin_flip,nBas,nC,nO,nV,nR,ENuc,EUH
ERI_aabb,ERI_bbbb,Om1sf,X1sf,Y1sf,&
Om2sf,X2sf,Y2sf,Ec_ppURPA(ispin))
call print_excitation('pp-RPA (N+2)',iblock,nP_sf,Om1sf)
call print_excitation('pp-RPA (N-2)',iblock,nH_sf,Om2sf)
call print_excitation_energies('ppRPA@HF (N+2)',iblock,nP_sf,Om1sf)
call print_excitation_energies('ppRPA@HF (N-2)',iblock,nH_sf,Om2sf)
write(*,*)
write(*,*)'-------------------------------------------------------------------------------'