PTEROSOR/QuAcK
10
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2024-11-05 05:33:50 +01:00
Code Releases Activity

ppRPA orlando orthogonalization

Browse Source
This commit is contained in:
Pierre-Francois Loos 2020-03-07 08:58:50 +01:00
parent c71ccda4df
commit 08d97376f2
15 changed files with 168 additions and 78 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
examples/basis.Hu
View File

@ -1,6 +1,6 @@
1 1 1 1
S 1 1.00 S 1 1.00
1.0000000 1.0000000 1 1.0000000 1.0000000
2 1 2 1
S 1 1.00 S 1 1.00
1.0000000 1.0000000 1 1.0000000 1.0000000

2
examples/molecule.F2
View File

@ -2,4 +2,4 @@
2 9 9 0 0 2 9 9 0 0
# Znuc x y z # Znuc x y z
F 0. 0. 0. F 0. 0. 0.
F 0. 0. 2.64 F 0. 0. 2.6682933

2
examples/molecule.LiH
View File

@ -2,4 +2,4 @@
2 2 2 0 0 2 2 2 0 0
# Znuc x y z # Znuc x y z
Li 0. 0. 0. Li 0. 0. 0.
H 0. 0. 2.994 H 0. 0. 3.0141132

67
input/basis
View File

@ -1,30 +1,49 @@
1 6 1 9
S 3 S 8
1 33.8700000 0.0060680 1 1469.0000000 0.0007660
2 5.0950000 0.0453080 2 220.5000000 0.0058920
3 1.1590000 0.2028220 3 50.2600000 0.0296710
4 14.2400000 0.1091800
5 4.5810000 0.2827890
6 1.5800000 0.4531230
7 0.5640000 0.2747740
8 0.0734500 0.0097510
S 8
1 1469.0000000 -0.0001200
2 220.5000000 -0.0009230
3 50.2600000 -0.0046890
4 14.2400000 -0.0176820
5 4.5810000 -0.0489020
6 1.5800000 -0.0960090
7 0.5640000 -0.1363800
8 0.0734500 0.5751020
S 1 S 1
1 0.3258000 1.0000000 1 0.0280500 1.0000000
S 1 S 1
1 0.1027000 1.0000000 1 0.0086400 1.0000000
P 3
1 1.5340000 0.0227840
2 0.2749000 0.1391070
3 0.0736200 0.5003750
P 1 P 1
1 1.4070000 1.0000000 1 0.0240300 1.0000000
P 1 P 1
1 0.3880000 1.0000000 1 0.0057900 1.0000000
D 1 D 1
1 1.0570000 1.0000000 1 0.1239000 1.0000000
2 6
S 3
1 33.8700000 0.0060680
2 5.0950000 0.0453080
3 1.1590000 0.2028220
S 1
1 0.3258000 1.0000000
S 1
1 0.1027000 1.0000000
P 1
1 1.4070000 1.0000000
P 1
1 0.3880000 1.0000000
D 1 D 1
1 1.0570000 1.0000000 1 0.0725000 1.0000000
2 5
S 3
1 13.0100000 0.0196850
2 1.9620000 0.1379770
3 0.4446000 0.4781480
S 1
1 0.1220000 1.0000000
S 1
1 0.0297400 1.0000000
P 1
1 0.7270000 1.0000000
P 1
1 0.1410000 1.0000000

6
input/methods
View File

@ -7,10 +7,10 @@
# CIS RPA RPAx ppRPA ADC # CIS RPA RPAx ppRPA ADC
F F F F F F F F F F
# GF2 GF3 # GF2 GF3
F F T F
# G0W0 evGW qsGW # G0W0 evGW qsGW
T F F
# G0T0 evGT qsGT
F F F F F F
# G0T0 evGT qsGT
T F F
# MCMP2 # MCMP2
F F

6
input/molecule
View File

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

4
input/molecule.xyz
View File

@ -1,4 +1,4 @@
2 2
H 0.0000000000 0.0000000000 0.0000000000 Li 0.0000000000 0.0000000000 0.0000000000
H 0.0000000000 0.0000000000 0.7408481486 H 0.0000000000 0.0000000000 1.5950001314

4
input/options
View File

@ -5,11 +5,11 @@
# CC: maxSCF thresh DIIS n_diis # CC: maxSCF thresh DIIS n_diis
64 0.0000001 F 1 64 0.0000001 F 1
# CIS/TDHF/BSE: singlet triplet # CIS/TDHF/BSE: singlet triplet
T F T T
# GF: maxSCF thresh DIIS n_diis renormalization # GF: maxSCF thresh DIIS n_diis renormalization
64 0.00001 T 5 3 64 0.00001 T 5 3
# GW: maxSCF thresh DIIS n_diis COHSEX SOSEX BSE TDA G0W GW0 lin eta # GW: maxSCF thresh DIIS n_diis COHSEX SOSEX BSE TDA G0W GW0 lin eta
256 0.00001 T 5 F F F F F F F 0.000 256 0.00001 T 5 F F T F F F F 0.000
# ACFDT: AC Kx XBS # ACFDT: AC Kx XBS
T F T T F T
# MCMP2: nMC nEq nWalk dt nPrint iSeed doDrift # MCMP2: nMC nEq nWalk dt nPrint iSeed doDrift

67
input/weight
View File

@ -1,30 +1,49 @@
1 6 1 9
S 3 S 8
1 33.8700000 0.0060680 1 1469.0000000 0.0007660
2 5.0950000 0.0453080 2 220.5000000 0.0058920
3 1.1590000 0.2028220 3 50.2600000 0.0296710
4 14.2400000 0.1091800
5 4.5810000 0.2827890
6 1.5800000 0.4531230
7 0.5640000 0.2747740
8 0.0734500 0.0097510
S 8
1 1469.0000000 -0.0001200
2 220.5000000 -0.0009230
3 50.2600000 -0.0046890
4 14.2400000 -0.0176820
5 4.5810000 -0.0489020
6 1.5800000 -0.0960090
7 0.5640000 -0.1363800
8 0.0734500 0.5751020
S 1 S 1
1 0.3258000 1.0000000 1 0.0280500 1.0000000
S 1 S 1
1 0.1027000 1.0000000 1 0.0086400 1.0000000
P 3
1 1.5340000 0.0227840
2 0.2749000 0.1391070
3 0.0736200 0.5003750
P 1 P 1
1 1.4070000 1.0000000 1 0.0240300 1.0000000
P 1 P 1
1 0.3880000 1.0000000 1 0.0057900 1.0000000
D 1 D 1
1 1.0570000 1.0000000 1 0.1239000 1.0000000
2 6
S 3
1 33.8700000 0.0060680
2 5.0950000 0.0453080
3 1.1590000 0.2028220
S 1
1 0.3258000 1.0000000
S 1
1 0.1027000 1.0000000
P 1
1 1.4070000 1.0000000
P 1
1 0.3880000 1.0000000
D 1 D 1
1 1.0570000 1.0000000 1 0.0725000 1.0000000
2 5
S 3
1 13.0100000 0.0196850
2 1.9620000 0.1379770
3 0.4446000 0.4781480
S 1
1 0.1220000 1.0000000
S 1
1 0.0297400 1.0000000
P 1
1 0.7270000 1.0000000
P 1
1 0.1410000 1.0000000

3
src/QuAcK/Bethe_Salpeter_A_matrix.f90
View File

@ -46,4 +46,7 @@ subroutine Bethe_Salpeter_A_matrix(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,Omega,rho,
enddo enddo
enddo enddo
print*,'BSE A'
call matout(nS,nS,A_lr)
end subroutine Bethe_Salpeter_A_matrix end subroutine Bethe_Salpeter_A_matrix

3
src/QuAcK/Bethe_Salpeter_B_matrix.f90
View File

@ -46,4 +46,7 @@ subroutine Bethe_Salpeter_B_matrix(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,Omega,rho,
enddo enddo
enddo enddo
print*,'BSE B'
call matout(nS,nS,B_lr)
end subroutine Bethe_Salpeter_B_matrix end subroutine Bethe_Salpeter_B_matrix

4
src/QuAcK/G0W0.f90
View File

@ -101,11 +101,11 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,COHSEX,SOSEX,BSE,TDA,singlet_manif
! Find all the roots of the QP equation if necessary ! Find all the roots of the QP equation if necessary
call QP_roots(nBas,nC,nO,nV,nR,nS,eta,eHF,Omega,rho,eGW) ! call QP_roots(nBas,nC,nO,nV,nR,nS,eta,eHF,Omega,rho,eGW)
! Dump results ! Dump results
! call print_excitation('RPA ',ispin,nS,Omega(:,ispin)) call print_excitation('RPA ',ispin,nS,Omega(:,ispin))
call print_G0W0(nBas,nO,eHF,ENuc,ERHF,SigC,Z,eGW,EcRPA(ispin),EcGM) call print_G0W0(nBas,nO,eHF,ENuc,ERHF,SigC,Z,eGW,EcRPA(ispin),EcGM)
! Compute the RPA correlation energy ! Compute the RPA correlation energy

6
src/QuAcK/linear_response_pp.f90
View File

@ -70,19 +70,19 @@ subroutine linear_response_pp(ispin,BSE,nBas,nC,nO,nV,nR,nOO,nVV,e,ERI,Omega1,X1
do ab=1,nVV do ab=1,nVV
do ij=1,nOO do ij=1,nOO
write(42,*) ab,ij,B(ab,ij) if(abs(B(ab,ij)) > 1d-15) write(42,*) ab,ij,B(ab,ij)
end do end do
end do end do
do ab=1,nVV do ab=1,nVV
do cd=1,nVV do cd=1,nVV
write(43,*) ab,cd,C(ab,cd) if(abs(C(ab,cd)) > 1d-15) write(43,*) ab,cd,C(ab,cd)
end do end do
end do end do
do ij=1,nOO do ij=1,nOO
do kl=1,nOO do kl=1,nOO
write(44,*) ij,kl,D(ij,kl) if(abs(D(ij,kl)) > 1d-15) write(44,*) ij,kl,D(ij,kl)
end do end do
end do end do

1
src/QuAcK/obj/.gitingore
View File

@ -1 +0,0 @@
*.o

55
src/QuAcK/sort_ppRPA.f90
View File

@ -15,6 +15,13 @@ subroutine sort_ppRPA(nOO,nVV,Omega,Z,Omega1,X1,Y1,Omega2,X2,Y2)
! Local variables ! Local variables
integer :: pq,ab,ij integer :: pq,ab,ij
double precision,allocatable :: M(:,:)
double precision,allocatable :: Z1(:,:)
double precision,allocatable :: Z2(:,:)
double precision,allocatable :: S1(:,:)
double precision,allocatable :: S2(:,:)
double precision,allocatable :: O1(:,:)
double precision,allocatable :: O2(:,:)
! Output variables ! Output variables
@ -25,6 +32,13 @@ subroutine sort_ppRPA(nOO,nVV,Omega,Z,Omega1,X1,Y1,Omega2,X2,Y2)
double precision,intent(out) :: X2(nVV,nOO) double precision,intent(out) :: X2(nVV,nOO)
double precision,intent(out) :: Y2(nOO,nOO) double precision,intent(out) :: Y2(nOO,nOO)
! Memory allocation
allocate(M(nOO+nVV,nOO+nVV), &
Z1(nOO+nVV,nVV),Z2(nOO+nVV,nOO), &
S1(nVV,nVV),S2(nOO,nOO), &
O1(nVV,nVV),O2(nOO,nOO))
! Initializatiom ! Initializatiom
Omega1(:) = 0d0 Omega1(:) = 0d0
@ -35,6 +49,20 @@ subroutine sort_ppRPA(nOO,nVV,Omega,Z,Omega1,X1,Y1,Omega2,X2,Y2)
X2(:,:) = 0d0 X2(:,:) = 0d0
Y2(:,:) = 0d0 Y2(:,:) = 0d0
! Compute metric
M(:,:) = 0d0
do ab=1,nVV
M(ab,ab) = 1d0
end do
do ij=1,nOO
M(nVV+ij,nVV+ij) = -1d0
end do
! Start sorting eigenvectors
ab = 0 ab = 0
ij = 0 ij = 0
@ -44,15 +72,17 @@ subroutine sort_ppRPA(nOO,nVV,Omega,Z,Omega1,X1,Y1,Omega2,X2,Y2)
ab = ab + 1 ab = ab + 1
Omega1(ab) = Omega(pq) Omega1(ab) = Omega(pq)
X1(1:nVV,ab) = Z( 1: nVV,pq) Z1(1:nOO+nVV,ab) = Z(1:nOO+nVV,pq)
Y1(1:nOO,ab) = Z(nVV+1:nOO+nVV,pq) ! X1(1:nVV,ab) = Z( 1: nVV,pq)
! Y1(1:nOO,ab) = Z(nVV+1:nOO+nVV,pq)
else else
ij = ij + 1 ij = ij + 1
Omega2(ij) = Omega(pq) Omega2(ij) = Omega(pq)
X2(1:nVV,ij) = Z( 1: nVV,pq) Z2(1:nOO+nVV,ij) = Z(1:nOO+nVV,pq)
Y2(1:nOO,ij) = Z(nVV+1:nOO+nVV,pq) ! X2(1:nVV,ij) = Z( 1: nVV,pq)
! Y2(1:nOO,ij) = Z(nVV+1:nOO+nVV,pq)
end if end if
@ -69,6 +99,23 @@ subroutine sort_ppRPA(nOO,nVV,Omega,Z,Omega1,X1,Y1,Omega2,X2,Y2)
call matout(nOO,1,Omega2(:)) call matout(nOO,1,Omega2(:))
write(*,*) write(*,*)
! Orthogonalize eigenvectors
S1 = matmul(transpose(Z1),matmul(M,Z1))
S2 = - matmul(transpose(Z2),matmul(M,Z2))
call orthogonalization_matrix(1,nVV,S1,O1)
call orthogonalization_matrix(1,nOO,S2,O2)
Z1 = matmul(Z1,O1)
Z2 = matmul(Z2,O2)
X1(1:nVV,1:nVV) = Z1( 1: nVV,1:nVV)
Y1(1:nOO,1:nVV) = Z1(nVV+1:nOO+nVV,1:nVV)
X2(1:nVV,1:nOO) = Z2( 1: nVV,1:nOO)
Y2(1:nOO,1:nOO) = Z2(nVV+1:nOO+nVV,1:nOO)
write(*,*) 'X1t.X1 - Y1t.Y1' write(*,*) 'X1t.X1 - Y1t.Y1'
call matout(nVV,nVV,matmul(transpose(X1),X1) - matmul(transpose(Y1),Y1)) call matout(nVV,nVV,matmul(transpose(X1),X1) - matmul(transpose(Y1),Y1))
write(*,*) 'X2t.X2 - Y2t.Y2' write(*,*) 'X2t.X2 - Y2t.Y2'