LCPQ/quack
4
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2024-12-23 04:43:53 +01:00
Code Releases Activity

fix DIIS in GW

Browse Source
This commit is contained in:
Pierre-Francois Loos 2019-04-01 15:53:39 +02:00
parent 89494730fa
commit 881220a97f
12 changed files with 197 additions and 69 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

100
input/basis
View File

@ -1,29 +1,73 @@
1 6 1 36
S 8 1.00
2940.0000000 0.0006800
441.2000000 0.0052360
100.5000000 0.0266060
28.4300000 0.0999930
9.1690000 0.2697020
3.1960000 0.4514690
1.1590000 0.2950740
0.1811000 0.0125870
S 8 1.00
2940.0000000 -0.0001230
441.2000000 -0.0009660
100.5000000 -0.0048310
28.4300000 -0.0193140
9.1690000 -0.0532800
3.1960000 -0.1207230
1.1590000 -0.1334350
0.1811000 0.5307670
S 1 1.00 S 1 1.00
0.0589000 1.0000000 1.0000000 1.0000000
P 3 1.00 S 1 1.00
3.6190000 0.0291110 1.0000000 1.0000000
0.7110000 0.1693650 S 1 1.00
0.1951000 0.5134580 1.0000000 1.0000000
P 1 1.00 S 1 1.00
0.0601800 1.0000000 1.0000000 1.0000000
D 1 1.00 S 1 1.00
0.2380000 1.0000000 1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000

2
input/methods
View File

@ -9,6 +9,6 @@
# GF2 GF3 # GF2 GF3
F F F F
# G0W0 evGW qsGW # G0W0 evGW qsGW
T T F T T T
# MCMP2 # MCMP2
F F

4
input/molecule
View File

@ -1,4 +1,4 @@
# nAt nEla nElb nCore nRyd # nAt nEla nElb nCore nRyd
1 2 2 0 0 1 1 1 0 0
# Znuc x y z # Znuc x y z
Be 0.0 0.0 0.0 X 0.0 0.0 0.0

60
input/weight
View File

@ -1,29 +1,39 @@
1 6 1 10
S 8 1.00 S 9 1.00
2940.0000000 0.0006800 6863.0000000 0.0002360
441.2000000 0.0052360 1030.0000000 0.0018260
100.5000000 0.0266060 234.7000000 0.0094520
28.4300000 0.0999930 66.5600000 0.0379570
9.1690000 0.2697020 21.6900000 0.1199650
3.1960000 0.4514690 7.7340000 0.2821620
1.1590000 0.2950740 2.9160000 0.4274040
0.1811000 0.0125870 1.1300000 0.2662780
S 8 1.00 0.1101000 -0.0072750
2940.0000000 -0.0001230 S 9 1.00
441.2000000 -0.0009660 6863.0000000 -0.0000430
100.5000000 -0.0048310 1030.0000000 -0.0003330
28.4300000 -0.0193140 234.7000000 -0.0017360
9.1690000 -0.0532800 66.5600000 -0.0070120
3.1960000 -0.1207230 21.6900000 -0.0231260
1.1590000 -0.1334350 7.7340000 -0.0581380
0.1811000 0.5307670 2.9160000 -0.1145560
1.1300000 -0.1359080
0.1101000 0.5774410
S 1 1.00 S 1 1.00
0.0589000 1.0000000 0.2577000 1.0000000
S 1 1.00
0.0440900 1.0000000
P 3 1.00 P 3 1.00
3.6190000 0.0291110 7.4360000 0.0107360
0.7110000 0.1693650 1.5770000 0.0628540
0.1951000 0.5134580 0.4352000 0.2481800
P 1 1.00 P 1 1.00
0.0601800 1.0000000 0.1438000 1.0000000
P 1 1.00
0.0499400 1.0000000
D 1 1.00 D 1 1.00
0.2380000 1.0000000 0.3480000 1.0000000
D 1 1.00
0.1803000 1.0000000
F 1 1.00
0.3250000 1.0000000

2
src/QuAcK/G0W0.f90
View File

@ -1,5 +1,5 @@
subroutine G0W0(COHSEX,SOSEX,BSE,TDA,singlet_manifold,triplet_manifold, & subroutine G0W0(COHSEX,SOSEX,BSE,TDA,singlet_manifold,triplet_manifold, &
nBas,nC,nO,nV,nR,nS,ENuc,ERHF,Hc,P,ERI_AO_basis,ERI_MO_basis,cHF,eHF,eG0W0) nBas,nC,nO,nV,nR,nS,ENuc,ERHF,Hc,ERI_AO_basis,ERI_MO_basis,P,cHF,eHF,eG0W0)
! Perform G0W0 calculation ! Perform G0W0 calculation

7
src/QuAcK/QuAcK.f90
View File

@ -405,7 +405,7 @@ program QuAcK
call cpu_time(start_G0W0) call cpu_time(start_G0W0)
call G0W0(COHSEX,SOSEX,BSE,TDA,singlet_manifold,triplet_manifold, & call G0W0(COHSEX,SOSEX,BSE,TDA,singlet_manifold,triplet_manifold, &
nBas,nC,nO,nV,nR,nS,ENuc,ERHF,Hc,PHF,ERI_AO_basis,ERI_MO_basis,cHF,eHF,eG0W0) nBas,nC,nO,nV,nR,nS,ENuc,ERHF,Hc,ERI_AO_basis,ERI_MO_basis,PHF,cHF,eHF,eG0W0)
call cpu_time(end_G0W0) call cpu_time(end_G0W0)
t_G0W0 = end_G0W0 - start_G0W0 t_G0W0 = end_G0W0 - start_G0W0
@ -421,9 +421,8 @@ program QuAcK
if(doevGW) then if(doevGW) then
call cpu_time(start_evGW) call cpu_time(start_evGW)
call evGW(maxSCF_GW,thresh_GW,n_diis_GW, & call evGW(maxSCF_GW,thresh_GW,n_diis_GW,COHSEX,SOSEX,BSE,TDA,G0W,GW0,singlet_manifold,triplet_manifold,linearize, &
COHSEX,SOSEX,BSE,TDA,G0W,GW0,singlet_manifold,triplet_manifold,linearize, & nBas,nC(1),nO(1),nV(1),nR(1),nS(1),ENuc,ERHF,Hc,ERI_AO_basis,ERI_MO_basis,PHF,cHF,eHF,eG0W0)
nBas,nC,nO,nV,nR,nS,ENuc,ERHF,Hc,ERI_AO_basis,ERI_MO_basis,PHF,cHF,eHF,eHF)
call cpu_time(end_evGW) call cpu_time(end_evGW)
t_evGW = end_evGW - start_evGW t_evGW = end_evGW - start_evGW

13
src/QuAcK/evGW.f90
View File

@ -20,7 +20,8 @@ subroutine evGW(maxSCF,thresh,max_diis,COHSEX,SOSEX,BSE,TDA,G0W,GW0,singlet_mani
logical :: dRPA,linear_mixing logical :: dRPA,linear_mixing
integer :: ispin,nSCF,n_diis integer :: ispin,nSCF,n_diis
double precision :: Conv,EcRPA,lambda double precision :: rcond
double precision :: Conv,EcRPA,EcGM,lambda
double precision,allocatable :: error_diis(:,:),e_diis(:,:) double precision,allocatable :: error_diis(:,:),e_diis(:,:)
double precision,allocatable :: eGW(:),eOld(:),Z(:) double precision,allocatable :: eGW(:),eOld(:),Z(:)
double precision,allocatable :: H(:,:),SigmaC(:) double precision,allocatable :: H(:,:),SigmaC(:)
@ -98,13 +99,13 @@ subroutine evGW(maxSCF,thresh,max_diis,COHSEX,SOSEX,BSE,TDA,G0W,GW0,singlet_mani
if(G0W) then if(G0W) then
call self_energy_correlation_diag(COHSEX,SOSEX,nBas,nC,nO,nV,nR,nS,eHF, & call self_energy_correlation_diag(COHSEX,SOSEX,nBas,nC,nO,nV,nR,nS,eHF, &
Omega(:,ispin),rho(:,:,:,ispin),rhox(:,:,:,ispin),SigmaC) Omega(:,ispin),rho(:,:,:,ispin),rhox(:,:,:,ispin),EcGM,SigmaC)
call renormalization_factor(SOSEX,nBas,nC,nO,nV,nR,nS,eHF,Omega(:,ispin),rho(:,:,:,ispin),rhox(:,:,:,ispin),Z) call renormalization_factor(SOSEX,nBas,nC,nO,nV,nR,nS,eHF,Omega(:,ispin),rho(:,:,:,ispin),rhox(:,:,:,ispin),Z)
else else
call self_energy_correlation_diag(COHSEX,SOSEX,nBas,nC,nO,nV,nR,nS,eGW, & call self_energy_correlation_diag(COHSEX,SOSEX,nBas,nC,nO,nV,nR,nS,eGW, &
Omega(:,ispin),rho(:,:,:,ispin),rhox(:,:,:,ispin),SigmaC) Omega(:,ispin),rho(:,:,:,ispin),rhox(:,:,:,ispin),EcGM,SigmaC)
call renormalization_factor(SOSEX,nBas,nC,nO,nV,nR,nS,eGW,Omega(:,ispin),rho(:,:,:,ispin),rhox(:,:,:,ispin),Z) call renormalization_factor(SOSEX,nBas,nC,nO,nV,nR,nS,eGW,Omega(:,ispin),rho(:,:,:,ispin),rhox(:,:,:,ispin),Z)
endif endif
@ -139,7 +140,11 @@ subroutine evGW(maxSCF,thresh,max_diis,COHSEX,SOSEX,BSE,TDA,G0W,GW0,singlet_mani
else else
n_diis = min(n_diis+1,max_diis) n_diis = min(n_diis+1,max_diis)
call DIIS_extrapolation(nBas,nBas,n_diis,error_diis,e_diis,eGW-eOld,eGW) call DIIS_extrapolation(rcond,nBas,nBas,n_diis,error_diis,e_diis,eGW-eOld,eGW)
! Reset DIIS if required
if(abs(rcond) < 1d-15) n_diis = 0
endif endif

65
src/QuAcK/excitation_density_RI.f90 Normal file
View File

@ -0,0 +1,65 @@
subroutine excitation_density_RI(nBas,nC,nO,nR,nS,c,G,XpY,rho)
! Compute excitation densities in the RI approximation
implicit none
! Input variables
integer,intent(in) :: nBas,nC,nO,nR,nS
double precision,intent(in) :: c(nBas,nBas),G(nBas,nBas,nBas,nBas),XpY(nS,nS)
! Local variables
double precision,allocatable :: scr(:,:,:)
integer :: mu,nu,la,si,ia,jb,x,y,j,b
! Output variables
double precision,intent(out) :: rho(nBas,nBas,nS)
! Memory allocation
allocate(scr(nBas,nBas,nS))
rho(:,:,:) = 0d0
do nu=1,nBas
do si=1,nBas
do ia=1,nS
jb = 0
do j=nC+1,nO
do b=nO+1,nBas-nR
jb = jb + 1
rho(nu,si,ia) = rho(nu,si,ia) + c(nu,j)*XpY(ia,jb)*c(si,b)
enddo
enddo
enddo
enddo
enddo
scr(:,:,:) = 0d0
do mu=1,nBas
do la=1,nBas
do ia=1,nS
do nu=1,nBas
do si=1,nBas
scr(mu,la,ia) = scr(mu,la,ia) + G(mu,nu,la,si)*rho(nu,si,ia)
enddo
enddo
enddo
enddo
enddo
rho(:,:,:) = 0d0
do ia=1,nS
do x=nC+1,nBas-nR
do y=nC+1,nBas-nR
do mu=1,nBas
do la=1,nBas
rho(x,y,ia) = rho(x,y,ia) + c(mu,x)*scr(mu,la,ia)*c(la,y)
enddo
enddo
enddo
enddo
enddo
end subroutine excitation_density_RI

7
src/QuAcK/qsGW.f90
View File

@ -22,6 +22,7 @@ subroutine qsGW(maxSCF,thresh,max_diis,COHSEX,SOSEX,BSE,TDA,G0W,GW0,singlet_mani
logical :: dRPA logical :: dRPA
integer :: nSCF,nBasSq,ispin,n_diis integer :: nSCF,nBasSq,ispin,n_diis
double precision :: EcRPA,EcGM,Conv double precision :: EcRPA,EcGM,Conv
double precision :: rcond
double precision,external :: trace_matrix double precision,external :: trace_matrix
double precision,allocatable :: error_diis(:,:),F_diis(:,:) double precision,allocatable :: error_diis(:,:),F_diis(:,:)
double precision,allocatable :: Omega(:,:),XpY(:,:,:),rho(:,:,:,:),rhox(:,:,:,:) double precision,allocatable :: Omega(:,:),XpY(:,:,:),rho(:,:,:,:),rhox(:,:,:,:)
@ -136,7 +137,11 @@ subroutine qsGW(maxSCF,thresh,max_diis,COHSEX,SOSEX,BSE,TDA,G0W,GW0,singlet_mani
! DIIS extrapolation ! DIIS extrapolation
n_diis = min(n_diis+1,max_diis) n_diis = min(n_diis+1,max_diis)
call DIIS_extrapolation(nBasSq,nBasSq,n_diis,error_diis,F_diis,error,F) call DIIS_extrapolation(rcond,nBasSq,nBasSq,n_diis,error_diis,F_diis,error,F)
! Reset DIIS if required
if(abs(rcond) < 1d-15) n_diis = 0
! Diagonalize Hamiltonian in AO basis ! Diagonalize Hamiltonian in AO basis

2
src/QuAcK/self_energy_correlation.f90
View File

@ -66,7 +66,7 @@ subroutine self_energy_correlation(COHSEX,SOSEX,nBas,nC,nO,nV,nR,nS,e,Omega,rho,
enddo enddo
enddo enddo
EcGM=0d0 EcGM = 0d0
do i=nC+1,nO do i=nC+1,nO
EcGM = EcGM + SigC(i,i) EcGM = EcGM + SigC(i,i)
enddo enddo

4
src/QuAcK/self_energy_correlation_diag.f90
View File

@ -66,7 +66,7 @@ subroutine self_energy_correlation_diag(COHSEX,SOSEX,nBas,nC,nO,nV,nR,nS,e,Omega
! GM correlation energy ! GM correlation energy
EcGM=0d0 EcGM = 0d0
do i=nC+1,nO do i=nC+1,nO
EcGM = EcGM + SigC(i) EcGM = EcGM + SigC(i)
enddo enddo
@ -109,7 +109,7 @@ subroutine self_energy_correlation_diag(COHSEX,SOSEX,nBas,nC,nO,nV,nR,nS,e,Omega
! GM correlation energy ! GM correlation energy
EcGM=0d0 EcGM = 0d0
do i=nC+1,nO do i=nC+1,nO
do a=nO+1,nBas-nR do a=nO+1,nBas-nR
jb = 0 jb = 0

0
src/QuAcK/DIIS_extrapolation.f90 → src/utils/DIIS_extrapolation.f90
View File