4
1
mirror of https://github.com/pfloos/quack synced 2024-06-25 22:52:18 +02:00

GW gap warning and antiresonnant term

This commit is contained in:
Pierre-Francois Loos 2020-05-15 13:08:07 +02:00
parent 431eab634d
commit 61832f1259
9 changed files with 148 additions and 123 deletions

View File

@ -1,62 +1,71 @@
1 14
S 3
1 82.6400000 0.0020060
2 12.4100000 0.0153430
3 2.8240000 0.0755790
1 9
S 8
1 9046.0000000 0.0007000
2 1357.0000000 0.0053890
3 309.3000000 0.0274060
4 87.7300000 0.1032070
5 28.5600000 0.2787230
6 10.2100000 0.4485400
7 3.8380000 0.2782380
8 0.7466000 0.0154400
S 8
1 9046.0000000 -0.0001530
2 1357.0000000 -0.0012080
3 309.3000000 -0.0059920
4 87.7300000 -0.0245440
5 28.5600000 -0.0674590
6 10.2100000 -0.1580780
7 3.8380000 -0.1218310
8 0.7466000 0.5490030
S 1
1 0.7977000 1.0000000
1 0.2248000 1.0000000
S 1
1 0.2581000 1.0000000
S 1
1 0.0898900 1.0000000
S 1
1 0.0236300 1.0000000
1 0.0612400 1.0000000
P 3
1 13.5500000 0.0399190
2 2.9170000 0.2171690
3 0.7973000 0.5103190
P 1
1 2.2920000 1.0000000
1 0.2185000 1.0000000
P 1
1 0.8380000 1.0000000
P 1
1 0.2920000 1.0000000
P 1
1 0.0848000 1.0000000
1 0.0561100 1.0000000
D 1
1 2.0620000 1.0000000
1 0.8170000 1.0000000
D 1
1 0.6620000 1.0000000
D 1
1 0.1900000 1.0000000
F 1
1 1.3970000 1.0000000
F 1
1 0.3600000 1.0000000
2 14
S 3
1 82.6400000 0.0020060
2 12.4100000 0.0153430
3 2.8240000 0.0755790
1 0.2300000 1.0000000
2 9
S 8
1 9046.0000000 0.0007000
2 1357.0000000 0.0053890
3 309.3000000 0.0274060
4 87.7300000 0.1032070
5 28.5600000 0.2787230
6 10.2100000 0.4485400
7 3.8380000 0.2782380
8 0.7466000 0.0154400
S 8
1 9046.0000000 -0.0001530
2 1357.0000000 -0.0012080
3 309.3000000 -0.0059920
4 87.7300000 -0.0245440
5 28.5600000 -0.0674590
6 10.2100000 -0.1580780
7 3.8380000 -0.1218310
8 0.7466000 0.5490030
S 1
1 0.7977000 1.0000000
1 0.2248000 1.0000000
S 1
1 0.2581000 1.0000000
S 1
1 0.0898900 1.0000000
S 1
1 0.0236300 1.0000000
1 0.0612400 1.0000000
P 3
1 13.5500000 0.0399190
2 2.9170000 0.2171690
3 0.7973000 0.5103190
P 1
1 2.2920000 1.0000000
1 0.2185000 1.0000000
P 1
1 0.8380000 1.0000000
P 1
1 0.2920000 1.0000000
P 1
1 0.0848000 1.0000000
1 0.0561100 1.0000000
D 1
1 2.0620000 1.0000000
1 0.8170000 1.0000000
D 1
1 0.6620000 1.0000000
D 1
1 0.1900000 1.0000000
F 1
1 1.3970000 1.0000000
F 1
1 0.3600000 1.0000000
1 0.2300000 1.0000000

View File

@ -13,7 +13,7 @@
# GGA = 2:
# Hybrid = 4:
# Hartree-Fock = 666
1 RMFL20
1 RVWN5
# quadrature grid SG-n
1
# Number of states in ensemble (nEns)

View File

@ -15,6 +15,6 @@
# G0W0 evGW qsGW
T F F
# G0T0 evGT qsGT
T F F
F F F
# MCMP2
F

View File

@ -1,5 +1,5 @@
# nAt nEla nElb nCore nRyd
2 1 1 0 0
2 7 7 0 0
# Znuc x y z
H 0. 0. 0.
H 0. 0. 1.4
N 0. 0. -1.04008632
N 0. 0. +1.04008632

View File

@ -1,4 +1,4 @@
2
H 0.0000000000 0.0000000000 0.0000000000
H 0.0000000000 0.0000000000 0.7408481486
N 0.0000000000 0.0000000000 -0.5503900175
N 0.0000000000 0.0000000000 0.5503900175

View File

@ -1,62 +1,71 @@
1 14
S 3
1 82.6400000 0.0020060
2 12.4100000 0.0153430
3 2.8240000 0.0755790
1 9
S 8
1 9046.0000000 0.0007000
2 1357.0000000 0.0053890
3 309.3000000 0.0274060
4 87.7300000 0.1032070
5 28.5600000 0.2787230
6 10.2100000 0.4485400
7 3.8380000 0.2782380
8 0.7466000 0.0154400
S 8
1 9046.0000000 -0.0001530
2 1357.0000000 -0.0012080
3 309.3000000 -0.0059920
4 87.7300000 -0.0245440
5 28.5600000 -0.0674590
6 10.2100000 -0.1580780
7 3.8380000 -0.1218310
8 0.7466000 0.5490030
S 1
1 0.7977000 1.0000000
1 0.2248000 1.0000000
S 1
1 0.2581000 1.0000000
S 1
1 0.0898900 1.0000000
S 1
1 0.0236300 1.0000000
1 0.0612400 1.0000000
P 3
1 13.5500000 0.0399190
2 2.9170000 0.2171690
3 0.7973000 0.5103190
P 1
1 2.2920000 1.0000000
1 0.2185000 1.0000000
P 1
1 0.8380000 1.0000000
P 1
1 0.2920000 1.0000000
P 1
1 0.0848000 1.0000000
1 0.0561100 1.0000000
D 1
1 2.0620000 1.0000000
1 0.8170000 1.0000000
D 1
1 0.6620000 1.0000000
D 1
1 0.1900000 1.0000000
F 1
1 1.3970000 1.0000000
F 1
1 0.3600000 1.0000000
2 14
S 3
1 82.6400000 0.0020060
2 12.4100000 0.0153430
3 2.8240000 0.0755790
1 0.2300000 1.0000000
2 9
S 8
1 9046.0000000 0.0007000
2 1357.0000000 0.0053890
3 309.3000000 0.0274060
4 87.7300000 0.1032070
5 28.5600000 0.2787230
6 10.2100000 0.4485400
7 3.8380000 0.2782380
8 0.7466000 0.0154400
S 8
1 9046.0000000 -0.0001530
2 1357.0000000 -0.0012080
3 309.3000000 -0.0059920
4 87.7300000 -0.0245440
5 28.5600000 -0.0674590
6 10.2100000 -0.1580780
7 3.8380000 -0.1218310
8 0.7466000 0.5490030
S 1
1 0.7977000 1.0000000
1 0.2248000 1.0000000
S 1
1 0.2581000 1.0000000
S 1
1 0.0898900 1.0000000
S 1
1 0.0236300 1.0000000
1 0.0612400 1.0000000
P 3
1 13.5500000 0.0399190
2 2.9170000 0.2171690
3 0.7973000 0.5103190
P 1
1 2.2920000 1.0000000
1 0.2185000 1.0000000
P 1
1 0.8380000 1.0000000
P 1
1 0.2920000 1.0000000
P 1
1 0.0848000 1.0000000
1 0.0561100 1.0000000
D 1
1 2.0620000 1.0000000
1 0.8170000 1.0000000
D 1
1 0.6620000 1.0000000
D 1
1 0.1900000 1.0000000
F 1
1 1.3970000 1.0000000
F 1
1 0.3600000 1.0000000
1 0.2300000 1.0000000

View File

@ -58,11 +58,11 @@ subroutine Bethe_Salpeter_B_matrix_dynamic(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,Om
chi = 0d0
do kc=1,maxS
eps = (OmBSE - OmRPA(kc) - (eGW(a) - eGW(i)))**2 + eta**2
chi = chi + rho(i,b,kc)*rho(a,j,kc)*(OmBSE - OmRPA(kc) - (eGW(a) - eGW(i)))/eps
eps = (OmBSE - OmRPA(kc) - (eGW(j) - eGW(i)))**2 + eta**2
chi = chi + rho(i,b,kc)*rho(a,j,kc)*(OmBSE - OmRPA(kc) - (eGW(j) - eGW(i)))/eps
eps = (OmBSE - OmRPA(kc) + (eGW(b) - eGW(j)))**2 + eta**2
chi = chi + rho(i,b,kc)*rho(a,j,kc)*(OmBSE - OmRPA(kc) + (eGW(b) - eGW(j)))/eps
eps = (OmBSE - OmRPA(kc) - (eGW(a) - eGW(b)))**2 + eta**2
chi = chi + rho(i,b,kc)*rho(a,j,kc)*(OmBSE - OmRPA(kc) - (eGW(a) - eGW(b)))/eps
enddo

View File

@ -27,6 +27,7 @@ subroutine Bethe_Salpeter_dynamic_perturbation(TDA,eta,nBas,nC,nO,nV,nR,nS,eGW,O
integer :: ia
integer,parameter :: maxS = 10
double precision :: gapGW
double precision,allocatable :: OmDyn(:)
double precision,allocatable :: ZDyn(:)
@ -38,30 +39,37 @@ subroutine Bethe_Salpeter_dynamic_perturbation(TDA,eta,nBas,nC,nO,nV,nR,nS,eGW,O
! Memory allocation
allocate(OmDyn(nS),ZDyn(nS),X(nS),Y(nS),A_dyn(nS,nS),Z_dyn(nS,nS))
allocate(OmDyn(nS),ZDyn(nS),X(nS),Y(nS),A_dyn(nS,nS),B_dyn(nS,nS),Z_dyn(nS,nS))
gapGW = eGW(nO+1) - eGW(nO)
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,*) ' First-order dynamical correction to static Bethe-Salpeter excitation energies '
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,'(2X,A5,1X,A20,1X,A20,1X,A20,1X,A20)') '#','Static (eV)','Dynamic (eV)','Correction (eV)','Renorm. (eV)'
write(*,*) '---------------------------------------------------------------------------------------------------'
do ia=1,min(nS,maxS)
X(:) = 0.5d0*(XpY(ia,:) + XmY(ia,:))
Y(:) = 0.5d0*(XpY(ia,:) - XmY(ia,:))
call Bethe_Salpeter_A_matrix_dynamic(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW(:),OmRPA(:),OmBSE(ia),rho(:,:,:),A_dyn(:,:))
! Renormalization factor
call Bethe_Salpeter_Z_matrix_dynamic(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW(:),OmRPA(:),OmBSE(ia),rho(:,:,:),Z_dyn(:,:))
ZDyn(ia) = dot_product(X(:),matmul(Z_dyn(:,:),X(:)))
ZDyn(ia) = 1d0/(1d0 - ZDyn(ia))
! First-order correction
if(.true.) then
ZDyn(ia) = dot_product(X(:),matmul(Z_dyn(:,:),X(:)))
OmDyn(ia) = dot_product(X(:),matmul(A_dyn(:,:),X(:)))
else
allocate(B_dyn(nS,nS))
call Bethe_Salpeter_B_matrix_dynamic(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW(:),OmRPA(:),OmBSE(ia),rho(:,:,:),B_dyn(:,:))
OmDyn(ia) = dot_product(X(:),matmul(A_dyn(:,:),X(:))) &
@ -71,14 +79,13 @@ subroutine Bethe_Salpeter_dynamic_perturbation(TDA,eta,nBas,nC,nO,nV,nR,nS,eGW,O
end if
! Renormalization factor
ZDyn(ia) = 1d0/(1d0 - ZDyn(ia))
OmDyn(ia) = ZDyn(ia)*dot_product(X(:),matmul(A_dyn(:,:),X(:)))
OmDyn(ia) = ZDyn(ia)*OmDyn(ia)
write(*,'(2X,I5,5X,F15.6,5X,F15.6,5X,F15.6,5X,F15.6)') &
ia,OmBSE(ia)*HaToeV,(OmBSE(ia)+OmDyn(ia))*HaToeV,OmDyn(ia)*HaToeV,ZDyn(ia)
if(OmBSE(ia) > gapGW) write(*,*) ' !!! BSE neutral excitation larger than the GW gap !!! '
end do
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,*)

View File

@ -82,8 +82,8 @@ subroutine LIM_RKS(x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,nGrid,weight, &
write(*,'(A40)') '*************************************************'
wLIM(1) = 0.5d0
wLIM(2) = 0.5d0
wLIM(3) = 0.0d0
wLIM(2) = 0.0d0
wLIM(3) = 0.5d0
do iEns=1,nEns
write(*,'(A20,I2,A2,F16.10)') ' Weight of state ',iEns,': ',wLIM(iEns)
@ -112,8 +112,8 @@ subroutine LIM_RKS(x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,nGrid,weight, &
write(*,'(A40)') '*************************************************'
write(*,*)
call GOK_RKS(.true.,x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,wLIM,nGrid,weight,maxSCF,thresh, &
max_diis,guess_type,nBas,AO,dAO,nO,nV,S,T,V,Hc,ERI,X,ENuc,Ew(3),c)
! call GOK_RKS(.true.,x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,wLIM,nGrid,weight,maxSCF,thresh, &
! max_diis,guess_type,nBas,AO,dAO,nO,nV,S,T,V,Hc,ERI,X,ENuc,Ew(3),c)
!------------------------------------------------------------------------
@ -121,7 +121,7 @@ subroutine LIM_RKS(x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,nGrid,weight, &
!------------------------------------------------------------------------
Om(2) = 2d0*(Ew(2) - Ew(1))
Om(3) = 3d0*(Ew(3) - Ew(2)) + 0.5d0*Om(2)
! Om(3) = 3d0*(Ew(3) - Ew(2)) + 0.5d0*Om(2)
write(*,'(A60)') '-------------------------------------------------'
write(*,'(A60)') ' LINEAR INTERPOLATION METHOD EXCITATION ENERGIES '