diff --git a/include/parameters.h b/include/parameters.h index d1d4fcb..d504995 100644 --- a/include/parameters.h +++ b/include/parameters.h @@ -1,6 +1,6 @@ integer,parameter :: ncart = 3 integer,parameter :: nspin = 2 - integer,parameter :: nsp = 3 + integer,parameter :: nsp = 3 integer,parameter :: maxEns = 10 integer,parameter :: maxShell = 512 integer,parameter :: maxL = 7 diff --git a/input/options b/input/options index 4dfc2da..5459e80 100644 --- a/input/options +++ b/input/options @@ -4,8 +4,8 @@ # CC: maxSCF thresh DIIS n_diis 64 0.0000001 T 5 -# spin: singlet triplet spin_conserved spin_flip TDA - T T T T F +# spin: TDA singlet triplet spin_conserved spin_flip + F T T T T # GF: maxSCF thresh DIIS n_diis lin eta renorm 256 0.00001 T 5 T 0.0 3 # GW/GT: maxSCF thresh DIIS n_diis lin eta COHSEX SOSEX TDA_W G0W GW0 @@ -13,6 +13,6 @@ # ACFDT: AC Kx XBS F F T # BSE: BSE dBSE dTDA evDyn - F F T F + T F T F # MCMP2: nMC nEq nWalk dt nPrint iSeed doDrift 1000000 100000 10 0.3 10000 1234 T diff --git a/src/CI/UCIS.f90 b/src/CI/UCIS.f90 index b846005..e9596e5 100644 --- a/src/CI/UCIS.f90 +++ b/src/CI/UCIS.f90 @@ -1,4 +1,4 @@ -subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab, & +subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,ERI_bbbb, & dipole_int_aa,dipole_int_bb,eHF,cHF,S) ! Perform configuration interaction single calculation` @@ -22,7 +22,6 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int_aa(nBas,nBas,ncart,nspin) double precision,intent(in) :: dipole_int_bb(nBas,nBas,ncart,nspin) @@ -72,7 +71,7 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E allocate(A_sc(nS_sc,nS_sc),Omega_sc(nS_sc)) call unrestricted_linear_response_A_matrix(ispin,.false.,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,eHF, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,A_sc) + ERI_aaaa,ERI_aabb,ERI_bbbb,A_sc) if(dump_matrix) then print*,'CIS matrix (spin-conserved transitions)' @@ -112,7 +111,7 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E allocate(A_sf(nS_sf,nS_sf),Omega_sf(nS_sf),Z_sf(nS_sf,nS_sf)) call unrestricted_linear_response_A_matrix(ispin,.false.,nBas,nC,nO,nV,nR,nS_ab,nS_ba,nS_sf,lambda,eHF, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,A_sf) + ERI_aaaa,ERI_aabb,ERI_bbbb,A_sf) if(dump_matrix) then print*,'CIS matrix (spin-conserved transitions)' diff --git a/src/utils/print_excitation.f90 b/src/LR/print_excitation.f90 similarity index 96% rename from src/utils/print_excitation.f90 rename to src/LR/print_excitation.f90 index c2bbdb4..26ee587 100644 --- a/src/utils/print_excitation.f90 +++ b/src/LR/print_excitation.f90 @@ -14,7 +14,7 @@ subroutine print_excitation(method,ispin,nS,Omega) ! Local variables character*14 :: spin_manifold - integer,parameter :: maxS = 32 + integer,parameter :: maxS = 10 integer :: ia if(ispin == 1) spin_manifold = 'singlet' diff --git a/src/LR/unrestricted_S2_expval.f90 b/src/LR/unrestricted_S2_expval.f90 index 559d6fc..6dc7018 100644 --- a/src/LR/unrestricted_S2_expval.f90 +++ b/src/LR/unrestricted_S2_expval.f90 @@ -152,23 +152,23 @@ subroutine unrestricted_S2_expval(ispin,nBas,nC,nO,nV,nR,nS,nSa,nSb,nSt,maxS,c,S Yat(:,:) = transpose(Ya(:,:)) Ybt(:,:) = transpose(Yb(:,:)) - S2(m) = 0d0 + S2(m) = S2_gs & -! S2(m) = S2_gs & -! + trace_matrix(nV(1),matmul(Xat,matmul(OO,matmul(OOt,Xa)))) & -! + trace_matrix(nV(2),matmul(Xbt,matmul(OOt,matmul(OO,Xb)))) & -! - trace_matrix(nO(1),matmul(Xa,matmul(VO,matmul(VOt,Xat)))) & -! - trace_matrix(nO(2),matmul(Xb,matmul(OVt,matmul(OV,Xbt)))) & -! - 2d0*trace_matrix(nO(1),matmul(OO,matmul(Xb,matmul(VVt,Xat)))) & -! -! - 2d0*trace_matrix(nV(2),matmul(OVt,matmul(Xa,matmul(VO,Yb)))) & -! - 2d0*trace_matrix(nV(1),matmul(VO,matmul(Xb,matmul(OVt,Ya)))) & -! -! - trace_matrix(nV(1),matmul(Yat,matmul(OO,matmul(OOt,Ya)))) & -! - trace_matrix(nV(2),matmul(Ybt,matmul(OOt,matmul(OO,Yb)))) & -! + trace_matrix(nO(1),matmul(Ya,matmul(VO,matmul(VOt,Yat)))) & -! + trace_matrix(nO(2),matmul(Yb,matmul(OVt,matmul(OV,Ybt)))) & -! + 2d0*trace_matrix(nO(1),matmul(Ya,matmul(VV,matmul(Ybt,OOt)))) + + trace_matrix(nV(1),matmul(Xbt,matmul(OOt,matmul(OO,Xb)))) & + - trace_matrix(nO(2),matmul(Xb,matmul(VO,matmul(VOt,Xbt)))) & + + trace_matrix(nO(2),matmul(Xb,VO))**2 & + + trace_matrix(nV(2),matmul(Yat,matmul(OO,matmul(OOt,Ya)))) & + + trace_matrix(nO(1),matmul(Ya,matmul(OVt,matmul(OV,Yat)))) & + + trace_matrix(nO(1),matmul(Ya,OVt))**2 & + - 2d0*trace_matrix(nO(2),matmul(Xb,VO))*trace_matrix(nO(1),matmul(Ya,OVt)) & + + + trace_matrix(nV(2),matmul(Xat,matmul(OO,matmul(OOt,Xa)))) & + - trace_matrix(nO(1),matmul(Xa,matmul(OVt,matmul(OV,Xat)))) & + + trace_matrix(nO(1),matmul(Xa,OVt))**2 & + + trace_matrix(nV(1),matmul(Ybt,matmul(OOt,matmul(OO,Yb)))) & + - trace_matrix(nO(2),matmul(Yb,matmul(VO,matmul(VOt,Ybt)))) & + + trace_matrix(nV(1),matmul(Ybt,VOt))**2 & + - 2d0*trace_matrix(nO(1),matmul(Xa,OVt))*trace_matrix(nO(2),matmul(Yb,VO)) end do diff --git a/src/LR/unrestricted_linear_response.f90 b/src/LR/unrestricted_linear_response.f90 index 3c747f8..32b2d6f 100644 --- a/src/LR/unrestricted_linear_response.f90 +++ b/src/LR/unrestricted_linear_response.f90 @@ -1,5 +1,5 @@ subroutine unrestricted_linear_response(ispin,dRPA,TDA,BSE,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_sc,lambda, & - e,ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,OmRPA,rho_RPA,EcRPA,Omega,XpY,XmY) + e,ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,EcRPA,Omega,XpY,XmY) ! Compute linear response for unrestricted formalism @@ -27,7 +27,6 @@ subroutine unrestricted_linear_response(ispin,dRPA,TDA,BSE,eta,nBas,nC,nO,nV,nR, double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) double precision,intent(in) :: OmRPA(nS_sc) double precision,intent(in) :: rho_RPA(nBas,nBas,nS_sc,nspin) @@ -58,11 +57,11 @@ subroutine unrestricted_linear_response(ispin,dRPA,TDA,BSE,eta,nBas,nC,nO,nV,nR, ! Build A and B matrices call unrestricted_linear_response_A_matrix(ispin,dRPA,nBas,nC,nO,nV,nR,nSa,nSb,nSt,lambda,e, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,A) + ERI_aaaa,ERI_aabb,ERI_bbbb,A) if(BSE) & call unrestricted_Bethe_Salpeter_A_matrix(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_sc,lambda,e, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,OmRPA,rho_RPA,A) + ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,A) ! Tamm-Dancoff approximation @@ -70,11 +69,11 @@ subroutine unrestricted_linear_response(ispin,dRPA,TDA,BSE,eta,nBas,nC,nO,nV,nR, if(.not. TDA) then call unrestricted_linear_response_B_matrix(ispin,dRPA,nBas,nC,nO,nV,nR,nSa,nSb,nSt,lambda, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,B) + ERI_aaaa,ERI_aabb,ERI_bbbb,B) if(BSE) & call unrestricted_Bethe_Salpeter_B_matrix(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_sc,lambda, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,OmRPA,rho_RPA,B) + ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,B) end if diff --git a/src/LR/unrestricted_linear_response_A_matrix.f90 b/src/LR/unrestricted_linear_response_A_matrix.f90 index 7878f9c..83cbf33 100644 --- a/src/LR/unrestricted_linear_response_A_matrix.f90 +++ b/src/LR/unrestricted_linear_response_A_matrix.f90 @@ -1,5 +1,5 @@ subroutine unrestricted_linear_response_A_matrix(ispin,dRPA,nBas,nC,nO,nV,nR,nSa,nSb,nSt,lambda, & - e,ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,A_lr) + e,ERI_aaaa,ERI_aabb,ERI_bbbb,A_lr) ! Compute linear response @@ -23,7 +23,6 @@ subroutine unrestricted_linear_response_A_matrix(ispin,dRPA,nBas,nC,nO,nV,nR,nSa double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) ! Local variables diff --git a/src/LR/unrestricted_linear_response_B_matrix.f90 b/src/LR/unrestricted_linear_response_B_matrix.f90 index d516743..fef8498 100644 --- a/src/LR/unrestricted_linear_response_B_matrix.f90 +++ b/src/LR/unrestricted_linear_response_B_matrix.f90 @@ -1,5 +1,5 @@ subroutine unrestricted_linear_response_B_matrix(ispin,dRPA,nBas,nC,nO,nV,nR,nSa,nSb,nSt,lambda, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,B_lr) + ERI_aaaa,ERI_aabb,ERI_bbbb,B_lr) ! Compute linear response @@ -22,7 +22,6 @@ subroutine unrestricted_linear_response_B_matrix(ispin,dRPA,nBas,nC,nO,nV,nR,nSa double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) ! Local variables diff --git a/src/MBPT/UG0W0.f90 b/src/MBPT/UG0W0.f90 index 53e3e41..843bf25 100644 --- a/src/MBPT/UG0W0.f90 +++ b/src/MBPT/UG0W0.f90 @@ -1,5 +1,5 @@ subroutine UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip, & - linearize,eta,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,Hc,ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab, & + linearize,eta,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,Hc,ERI_aaaa,ERI_aabb,ERI_bbbb, & dipole_int_aa,dipole_int_bb,PHF,cHF,eHF,eGW) ! Perform unrestricted G0W0 calculation @@ -40,7 +40,6 @@ subroutine UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,ev double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int_aa(nBas,nBas,ncart) double precision,intent(in) :: dipole_int_bb(nBas,nBas,ncart) @@ -118,7 +117,7 @@ subroutine UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,ev ispin = 1 call unrestricted_linear_response(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,ERI_abab,OmRPA,rho_RPA,EcRPA,OmRPA,XpY_RPA,XmY_RPA) + 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) @@ -167,7 +166,7 @@ subroutine UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,ev ! Compute RPA correlation energy call unrestricted_linear_response(ispin,.true.,TDA_W,.false.,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0, & - eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,OmRPA,rho_RPA,EcRPA,OmRPA,XpY_RPA,XmY_RPA) + eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,EcRPA,OmRPA,XpY_RPA,XmY_RPA) ! Dump results @@ -182,7 +181,7 @@ subroutine UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,ev if(BSE) then call unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta,nBas,nC,nO,nV,nR,nS, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,dipole_int_aa,dipole_int_bb,eHF,eGW,EcBSE) + ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,eHF,eGW,EcBSE) ! if(exchange_kernel) then ! diff --git a/src/MBPT/evUGW.f90 b/src/MBPT/evUGW.f90 index 322e13e..3323394 100644 --- a/src/MBPT/evUGW.f90 +++ b/src/MBPT/evUGW.f90 @@ -1,6 +1,6 @@ subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA, & G0W,GW0,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta,nBas,nC,nO,nV,nR,nS,ENuc, & - ERHF,Hc,ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF,eG0W0) + ERHF,Hc,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF,eG0W0) ! Perform self-consistent eigenvalue-only GW calculation @@ -45,7 +45,6 @@ subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int_aa(nBas,nBas,ncart) double precision,intent(in) :: dipole_int_bb(nBas,nBas,ncart) @@ -155,7 +154,7 @@ subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE if(.not. GW0 .or. nSCF == 0) then call unrestricted_linear_response(ispin,.true.,TDA_W,.false.,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0, & - eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,OmRPA,rho_RPA,EcRPA,OmRPA,XpY_RPA,XmY_RPA) + eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,EcRPA,OmRPA,XpY_RPA,XmY_RPA) endif @@ -256,7 +255,7 @@ subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE if(BSE) then call unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta,nBas,nC,nO,nV,nR,nS, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,dipole_int_aa,dipole_int_bb,eGW,eGW,EcBSE) + ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,eGW,eGW,EcBSE) ! if(exchange_kernel) then diff --git a/src/MBPT/unrestricted_Bethe_Salpeter.f90 b/src/MBPT/unrestricted_Bethe_Salpeter.f90 index d61182e..ac8bf1c 100644 --- a/src/MBPT/unrestricted_Bethe_Salpeter.f90 +++ b/src/MBPT/unrestricted_Bethe_Salpeter.f90 @@ -1,5 +1,5 @@ subroutine unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta, & - nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab, & + nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,ERI_bbbb, & dipole_int_aa,dipole_int_bb,eW,eGW,EcBSE) ! Compute the Bethe-Salpeter excitation energies @@ -29,7 +29,6 @@ subroutine unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved, double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int_aa(nBas,nBas,ncart) double precision,intent(in) :: dipole_int_bb(nBas,nBas,ncart) @@ -76,7 +75,7 @@ subroutine unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved, ! Compute spin-conserved RPA screening call unrestricted_linear_response(isp_W,.true.,TDA_W,.false.,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0, & - eW,ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,OmRPA,rho_RPA,EcRPA,OmRPA,XpY_RPA,XmY_RPA) + eW,ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,EcRPA,OmRPA,XpY_RPA,XmY_RPA) call unrestricted_excitation_density(nBas,nC,nO,nR,nS_aa,nS_bb,nS_sc,ERI_aaaa,ERI_aabb,ERI_bbbb,XpY_RPA,rho_RPA) @@ -95,7 +94,7 @@ subroutine unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved, ! Compute spin-conserved BSE excitation energies call unrestricted_linear_response(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,ERI_abab,OmRPA,rho_RPA,EcBSE(ispin), & + eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,EcBSE(ispin), & OmBSE_sc,XpY_BSE_sc,XmY_BSE_sc) call print_excitation('BSE@UG0W0',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, & @@ -107,7 +106,7 @@ subroutine unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved, if(dBSE) & call unrestricted_Bethe_Salpeter_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nS_aa,nS_bb,nS_sc,nS_sc, & - eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,dipole_int_aa,dipole_int_bb, & + eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb, & OmRPA,rho_RPA,OmBSE_sc,XpY_BSE_sc,XmY_BSE_sc) deallocate(OmBSE_sc,XpY_BSE_sc,XmY_BSE_sc) @@ -134,7 +133,7 @@ subroutine unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved, ! Compute spin-flip BSE excitation energies call unrestricted_linear_response(ispin,.true.,TDA,.true.,eta,nBas,nC,nO,nV,nR,nS_ab,nS_ba,nS_sf,nS_sc,1d0, & - eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,OmRPA,rho_RPA,EcBSE(ispin), & + eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,EcBSE(ispin), & OmBSE_sf,XpY_BSE_sf,XmY_BSE_sf) call print_excitation('BSE@UG0W0',6,nS_sf,OmBSE_sf) @@ -145,7 +144,7 @@ subroutine unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved, if(dBSE) & call unrestricted_Bethe_Salpeter_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nS_ab,nS_ba,nS_sf,nS_sc, & - eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,dipole_int_aa,dipole_int_bb, & + eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb, & OmRPA,rho_RPA,OmBSE_sf,XpY_BSE_sf,XmY_BSE_sf) deallocate(OmBSE_sf,XpY_BSE_sf,XmY_BSE_sf) diff --git a/src/MBPT/unrestricted_Bethe_Salpeter_A_matrix.f90 b/src/MBPT/unrestricted_Bethe_Salpeter_A_matrix.f90 index 31a6ef3..48752f8 100644 --- a/src/MBPT/unrestricted_Bethe_Salpeter_A_matrix.f90 +++ b/src/MBPT/unrestricted_Bethe_Salpeter_A_matrix.f90 @@ -1,5 +1,5 @@ subroutine unrestricted_Bethe_Salpeter_A_matrix(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_sc,lambda,eGW, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,Omega,rho,A_lr) + ERI_aaaa,ERI_aabb,ERI_bbbb,Omega,rho,A_lr) ! Compute the extra term for Bethe-Salpeter equation for linear response in the unrestricted formalism @@ -24,7 +24,6 @@ subroutine unrestricted_Bethe_Salpeter_A_matrix(ispin,eta,nBas,nC,nO,nV,nR,nSa,n double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) double precision,intent(in) :: Omega(nS_sc) double precision,intent(in) :: rho(nBas,nBas,nS_sc,nspin) diff --git a/src/MBPT/unrestricted_Bethe_Salpeter_A_matrix_dynamic.f90 b/src/MBPT/unrestricted_Bethe_Salpeter_A_matrix_dynamic.f90 index b985f23..c8e502c 100644 --- a/src/MBPT/unrestricted_Bethe_Salpeter_A_matrix_dynamic.f90 +++ b/src/MBPT/unrestricted_Bethe_Salpeter_A_matrix_dynamic.f90 @@ -1,5 +1,5 @@ subroutine unrestricted_Bethe_Salpeter_A_matrix_dynamic(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_sc,lambda,eGW, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,OmRPA,rho_RPA,OmBSE,A_dyn) + ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,OmBSE,A_dyn) ! Compute the extra term for dynamical Bethe-Salpeter equation for linear response in the unrestricted formalism @@ -24,7 +24,6 @@ subroutine unrestricted_Bethe_Salpeter_A_matrix_dynamic(ispin,eta,nBas,nC,nO,nV, double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) double precision,intent(in) :: OmRPA(nS_sc) double precision,intent(in) :: rho_RPA(nBas,nBas,nS_sc,nspin) double precision,intent(in) :: OmBSE diff --git a/src/MBPT/unrestricted_Bethe_Salpeter_B_matrix.f90 b/src/MBPT/unrestricted_Bethe_Salpeter_B_matrix.f90 index 69965f2..c4a45f5 100644 --- a/src/MBPT/unrestricted_Bethe_Salpeter_B_matrix.f90 +++ b/src/MBPT/unrestricted_Bethe_Salpeter_B_matrix.f90 @@ -1,5 +1,5 @@ subroutine unrestricted_Bethe_Salpeter_B_matrix(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_sc,lambda, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,Omega,rho,B_lr) + ERI_aaaa,ERI_aabb,ERI_bbbb,Omega,rho,B_lr) ! Compute the extra term for Bethe-Salpeter equation for linear response @@ -23,7 +23,6 @@ subroutine unrestricted_Bethe_Salpeter_B_matrix(ispin,eta,nBas,nC,nO,nV,nR,nSa,n double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) double precision,intent(in) :: Omega(nS_sc) double precision,intent(in) :: rho(nBas,nBas,nS_sc,nspin) diff --git a/src/MBPT/unrestricted_Bethe_Salpeter_ZA_matrix_dynamic.f90 b/src/MBPT/unrestricted_Bethe_Salpeter_ZA_matrix_dynamic.f90 index a811843..be1da42 100644 --- a/src/MBPT/unrestricted_Bethe_Salpeter_ZA_matrix_dynamic.f90 +++ b/src/MBPT/unrestricted_Bethe_Salpeter_ZA_matrix_dynamic.f90 @@ -1,5 +1,5 @@ subroutine unrestricted_Bethe_Salpeter_ZA_matrix_dynamic(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_sc,lambda,eGW, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,OmRPA,rho_RPA,OmBSE,ZA_dyn) + ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,OmBSE,ZA_dyn) ! Compute the extra term for dynamical Bethe-Salpeter equation for linear response in the unrestricted formalism @@ -24,7 +24,6 @@ subroutine unrestricted_Bethe_Salpeter_ZA_matrix_dynamic(ispin,eta,nBas,nC,nO,nV double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) double precision,intent(in) :: OmRPA(nS_sc) double precision,intent(in) :: rho_RPA(nBas,nBas,nS_sc,nspin) double precision,intent(in) :: OmBSE diff --git a/src/MBPT/unrestricted_Bethe_Salpeter_dynamic_perturbation.f90 b/src/MBPT/unrestricted_Bethe_Salpeter_dynamic_perturbation.f90 index 1057eee..03660ee 100644 --- a/src/MBPT/unrestricted_Bethe_Salpeter_dynamic_perturbation.f90 +++ b/src/MBPT/unrestricted_Bethe_Salpeter_dynamic_perturbation.f90 @@ -1,5 +1,5 @@ subroutine unrestricted_Bethe_Salpeter_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nSa,nSb,nSt,nS_sc,eGW, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,dipole_int_aa,dipole_int_bb, & + ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb, & OmRPA,rho_RPA,OmBSE,XpY_BSE,XmY_BSE) ! Compute dynamical effects via perturbation theory for BSE @@ -27,7 +27,6 @@ subroutine unrestricted_Bethe_Salpeter_dynamic_perturbation(ispin,dTDA,eta,nBas, double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int_aa(nBas,nBas,ncart) double precision,intent(in) :: dipole_int_bb(nBas,nBas,ncart) double precision,intent(in) :: OmRPA(nS_sc) @@ -85,12 +84,12 @@ subroutine unrestricted_Bethe_Salpeter_dynamic_perturbation(ispin,dTDA,eta,nBas, ! Resonant part of the BSE correction for dynamical TDA call unrestricted_Bethe_Salpeter_A_matrix_dynamic(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_sc,1d0,eGW, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,OmRPA,rho_RPA,OmBSE(ia),A_dyn) + ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,OmBSE(ia),A_dyn) ! Renormalization factor of the resonant parts for dynamical TDA call unrestricted_Bethe_Salpeter_ZA_matrix_dynamic(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_sc,1d0,eGW, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,OmRPA,rho_RPA,OmBSE(ia),ZA_dyn) + ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,OmBSE(ia),ZA_dyn) ZDyn(ia) = dot_product(X,matmul(ZA_dyn,X)) OmDyn(ia) = dot_product(X,matmul( A_dyn,X)) diff --git a/src/QuAcK/QuAcK.f90 b/src/QuAcK/QuAcK.f90 index 27206bf..58de284 100644 --- a/src/QuAcK/QuAcK.f90 +++ b/src/QuAcK/QuAcK.f90 @@ -70,7 +70,6 @@ program QuAcK double precision,allocatable :: ERI_MO_aaaa(:,:,:,:) double precision,allocatable :: ERI_MO_aabb(:,:,:,:) double precision,allocatable :: ERI_MO_bbbb(:,:,:,:) - double precision,allocatable :: ERI_MO_abab(:,:,:,:) double precision,allocatable :: ERI_ERF_AO(:,:,:,:) double precision,allocatable :: ERI_ERF_MO(:,:,:,:) double precision,allocatable :: F12(:,:,:,:),Yuk(:,:,:,:),FC(:,:,:,:,:,:) @@ -169,7 +168,7 @@ program QuAcK call read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_type, & maxSCF_CC,thresh_CC,DIIS_CC,n_diis_CC, & - singlet,triplet,spin_conserved,spin_flip,TDA, & + TDA,singlet,triplet,spin_conserved,spin_flip, & maxSCF_GF,thresh_GF,DIIS_GF,n_diis_GF,linGF,eta_GF,renormGF, & maxSCF_GW,thresh_GW,DIIS_GW,n_diis_GW,linGW,eta_GW, & COHSEX,SOSEX,TDA_W,G0W,GW0, & @@ -356,8 +355,7 @@ program QuAcK ! Memory allocation - allocate(ERI_MO_aaaa(nBas,nBas,nBas,nBas),ERI_MO_aabb(nBas,nBas,nBas,nBas), & - ERI_MO_bbbb(nBas,nBas,nBas,nBas),ERI_MO_abab(nBas,nBas,nBas,nBas)) + allocate(ERI_MO_aaaa(nBas,nBas,nBas,nBas),ERI_MO_aabb(nBas,nBas,nBas,nBas),ERI_MO_bbbb(nBas,nBas,nBas,nBas)) ! 4-index transform for (aa|aa) block @@ -383,14 +381,6 @@ program QuAcK ket2 = 2 call AOtoMO_integral_transform(bra1,bra2,ket1,ket2,nBas,cHF,ERI_AO,ERI_MO_bbbb) - ! 4-index transform for (ab|ab) block - - bra1 = 1 - bra2 = 2 - ket1 = 1 - ket2 = 2 - call AOtoMO_integral_transform(bra1,bra2,ket1,ket2,nBas,cHF,ERI_AO,ERI_MO_abab) - else ! Memory allocation @@ -612,7 +602,7 @@ program QuAcK if(unrestricted) then call UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_MO_aaaa,ERI_MO_aabb, & - ERI_MO_bbbb,ERI_MO_abab,dipole_int_aa,dipole_int_bb,eHF,cHF,S) + ERI_MO_bbbb,dipole_int_aa,dipole_int_bb,eHF,cHF,S) else @@ -669,7 +659,7 @@ program QuAcK if(unrestricted) then call UdRPA(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,0d0,nBas,nC,nO,nV,nR,nS,ENuc,EUHF, & - ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,ERI_MO_abab,dipole_int_aa,dipole_int_bb,eHF,cHF,S) + ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,dipole_int_aa,dipole_int_bb,eHF,cHF,S) else @@ -694,7 +684,7 @@ program QuAcK if(unrestricted) then call URPAx(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,0d0,nBas,nC,nO,nV,nR,nS,ENuc,EUHF, & - ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,ERI_MO_abab,dipole_int_aa,dipole_int_bb,eHF,cHF,S) + ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,dipole_int_aa,dipole_int_bb,eHF,cHF,S) else @@ -821,7 +811,7 @@ program QuAcK if(unrestricted) then call UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip, & - linGW,eta_GW,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,Hc,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,ERI_MO_abab, & + linGW,eta_GW,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,Hc,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb, & dipole_int_aa,dipole_int_bb,PHF,cHF,eHF,eG0W0) else @@ -849,7 +839,7 @@ program QuAcK call evUGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA, & G0W,GW0,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta_GW,nBas,nC,nO,nV,nR,nS,ENuc, & - EUHF,Hc,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,ERI_MO_abab,dipole_int_aa,dipole_int_bb, & + EUHF,Hc,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,dipole_int_aa,dipole_int_bb, & PHF,cHF,eHF,eG0W0) else diff --git a/src/QuAcK/read_options.f90 b/src/QuAcK/read_options.f90 index c9e72f3..a667a3c 100644 --- a/src/QuAcK/read_options.f90 +++ b/src/QuAcK/read_options.f90 @@ -1,11 +1,11 @@ -subroutine read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_type, & - maxSCF_CC,thresh_CC,DIIS_CC,n_diis_CC, & - singlet,triplet,spin_conserved,spin_flip,TDA, & - maxSCF_GF,thresh_GF,DIIS_GF,n_diis_GF,linGF,eta_GF,renormGF, & - maxSCF_GW,thresh_GW,DIIS_GW,n_diis_GW,linGW,eta_GW, & - COHSEX,SOSEX,TDA_W,G0W,GW0, & - doACFDT,exchange_kernel,doXBS, & - BSE,dBSE,dTDA,evDyn, & +subroutine read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_type, & + maxSCF_CC,thresh_CC,DIIS_CC,n_diis_CC, & + TDA,singlet,triplet,spin_conserved,spin_flip, & + maxSCF_GF,thresh_GF,DIIS_GF,n_diis_GF,linGF,eta_GF,renormGF, & + maxSCF_GW,thresh_GW,DIIS_GW,n_diis_GW,linGW,eta_GW, & + COHSEX,SOSEX,TDA_W,G0W,GW0, & + doACFDT,exchange_kernel,doXBS, & + BSE,dBSE,dTDA,evDyn, & nMC,nEq,nWalk,dt,nPrint,iSeed,doDrift) ! Read desired methods @@ -26,11 +26,11 @@ subroutine read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_t logical,intent(out) :: DIIS_CC integer,intent(out) :: n_diis_CC + logical,intent(out) :: TDA logical,intent(out) :: singlet logical,intent(out) :: triplet logical,intent(out) :: spin_conserved logical,intent(out) :: spin_flip - logical,intent(out) :: TDA integer,intent(out) :: maxSCF_GF double precision,intent(out) :: thresh_GF @@ -115,20 +115,20 @@ subroutine read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_t ! Read excited state options + TDA = .false. singlet = .false. triplet = .false. spin_conserved = .false. spin_flip = .false. - TDA = .false. read(1,*) read(1,*) answer1,answer2,answer3,answer4,answer5 - if(answer1 == 'T') singlet = .true. - if(answer2 == 'T') triplet = .true. - if(answer3 == 'T') spin_conserved = .true. - if(answer4 == 'T') spin_flip = .true. - if(answer5 == 'T') TDA = .true. + if(answer1 == 'T') TDA = .true. + if(answer2 == 'T') singlet = .true. + if(answer3 == 'T') triplet = .true. + if(answer4 == 'T') spin_conserved = .true. + if(answer5 == 'T') spin_flip = .true. ! Read Green function options diff --git a/src/RPA/URPAx.f90 b/src/RPA/URPAx.f90 index 33da34a..42cc633 100644 --- a/src/RPA/URPAx.f90 +++ b/src/RPA/URPAx.f90 @@ -1,5 +1,5 @@ subroutine URPAx(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,eta,nBas,nC,nO,nV,nR,nS,ENuc,EUHF, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,dipole_int_aa,dipole_int_bb,e,c,S) + ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,e,c,S) ! Perform random phase approximation calculation with exchange (aka TDHF) in the unrestricted formalism @@ -29,7 +29,6 @@ subroutine URPAx(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,eta,nBas,n double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int_aa(nBas,nBas,ncart) double precision,intent(in) :: dipole_int_bb(nBas,nBas,ncart) @@ -87,7 +86,7 @@ subroutine URPAx(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,eta,nBas,n allocate(Omega_sc(nS_sc),XpY_sc(nS_sc,nS_sc),XmY_sc(nS_sc,nS_sc)) call unrestricted_linear_response(ispin,.false.,TDA,.false.,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0,e, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,Omega_sc,rho_sc,EcRPAx(ispin),Omega_sc,XpY_sc,XmY_sc) + ERI_aaaa,ERI_aabb,ERI_bbbb,Omega_sc,rho_sc,EcRPAx(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) @@ -111,7 +110,7 @@ subroutine URPAx(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,eta,nBas,n allocate(Omega_sf(nS_sf),XpY_sf(nS_sf,nS_sf),XmY_sf(nS_sf,nS_sf)) call unrestricted_linear_response(ispin,.false.,TDA,.false.,eta,nBas,nC,nO,nV,nR,nS_ab,nS_ba,nS_sf,nS_sf,1d0,e, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,Omega_sf,rho_sf,EcRPAx(ispin),Omega_sf,XpY_sf,XmY_sf) + ERI_aaaa,ERI_aabb,ERI_bbbb,Omega_sf,rho_sf,EcRPAx(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) diff --git a/src/RPA/UdRPA.f90 b/src/RPA/UdRPA.f90 index 1becf23..22f29ec 100644 --- a/src/RPA/UdRPA.f90 +++ b/src/RPA/UdRPA.f90 @@ -1,5 +1,5 @@ subroutine UdRPA(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,eta,nBas,nC,nO,nV,nR,nS,ENuc,EUHF, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,dipole_int_aa,dipole_int_bb,e,c,S) + ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,e,c,S) ! Perform random phase approximation calculation with exchange (aka TDHF) in the unrestricted formalism @@ -29,7 +29,6 @@ subroutine UdRPA(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,eta,nBas,n double precision,intent(in) :: ERI_aaaa(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_aabb(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI_bbbb(nBas,nBas,nBas,nBas) - double precision,intent(in) :: ERI_abab(nBas,nBas,nBas,nBas) double precision,intent(in) :: dipole_int_aa(nBas,nBas,ncart) double precision,intent(in) :: dipole_int_bb(nBas,nBas,ncart) @@ -86,7 +85,7 @@ subroutine UdRPA(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,eta,nBas,n allocate(Omega_sc(nS_sc),XpY_sc(nS_sc,nS_sc),XmY_sc(nS_sc,nS_sc)) call unrestricted_linear_response(ispin,.true.,TDA,.false.,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0,e, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,Omega_sc,rho_sc,EcRPA(ispin),Omega_sc,XpY_sc,XmY_sc) + 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) @@ -109,7 +108,7 @@ subroutine UdRPA(TDA,doACFDT,exchange_kernel,spin_conserved,spin_flip,eta,nBas,n allocate(Omega_sf(nS_sf),XpY_sf(nS_sf,nS_sf),XmY_sf(nS_sf,nS_sf)) call unrestricted_linear_response(ispin,.true.,TDA,.false.,eta,nBas,nC,nO,nV,nR,nS_ab,nS_ba,nS_sf,nS_sf,1d0,e, & - ERI_aaaa,ERI_aabb,ERI_bbbb,ERI_abab,Omega_sf,rho_sf,EcRPA(ispin),Omega_sf,XpY_sf,XmY_sf) + 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)