10
1
mirror of https://github.com/pfloos/quack synced 2024-11-04 21:23:55 +01:00

cleaning ppRPA

This commit is contained in:
Pierre-Francois Loos 2020-04-15 08:53:02 +02:00
parent cfc910fb26
commit b5e4ded668
7 changed files with 71 additions and 65 deletions

View File

@ -1,9 +1,35 @@
1 3 1 9
S 3 S 8
1 38.3600000 0.0238090 1 2940.0000000 0.0006800
2 5.7700000 0.1548910 2 441.2000000 0.0052360
3 1.2400000 0.4699870 3 100.5000000 0.0266060
4 28.4300000 0.0999930
5 9.1690000 0.2697020
6 3.1960000 0.4514690
7 1.1590000 0.2950740
8 0.1811000 0.0125870
S 8
1 2940.0000000 -0.0001230
2 441.2000000 -0.0009660
3 100.5000000 -0.0048310
4 28.4300000 -0.0193140
5 9.1690000 -0.0532800
6 3.1960000 -0.1207230
7 1.1590000 -0.1334350
8 0.1811000 0.5307670
S 1 S 1
1 0.2976000 1.0000000 1 0.0589000 1.0000000
S 1
1 0.0187700 1.0000000
P 3
1 3.6190000 0.0291110
2 0.7110000 0.1693650
3 0.1951000 0.5134580
P 1 P 1
1 1.2750000 1.0000000 1 0.0601800 1.0000000
P 1
1 0.0085000 1.0000000
D 1
1 0.2380000 1.0000000
D 1
1 0.0740000 1.0000000

View File

@ -7,12 +7,12 @@
# drCCD rCCD lCCD pCCD # drCCD rCCD lCCD pCCD
F F F F F F F F
# CIS RPA RPAx ppRPA ADC # CIS RPA RPAx ppRPA ADC
F F F F F F F F T F
# G0F2 evGF2 G0F3 evGF3 # G0F2 evGF2 G0F3 evGF3
F F F F F F F F
# G0W0 evGW qsGW # G0W0 evGW qsGW
F F F F F F
# G0T0 evGT qsGT # G0T0 evGT qsGT
T F F F F F
# MCMP2 # MCMP2
F F

View File

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

View File

@ -1,3 +1,3 @@
1 1
He 0.0000000000 0.0000000000 0.0000000000 Be 0.0000000000 0.0000000000 0.0000000000

View File

@ -1,9 +1,35 @@
1 3 1 9
S 3 S 8
1 38.3600000 0.0238090 1 2940.0000000 0.0006800
2 5.7700000 0.1548910 2 441.2000000 0.0052360
3 1.2400000 0.4699870 3 100.5000000 0.0266060
4 28.4300000 0.0999930
5 9.1690000 0.2697020
6 3.1960000 0.4514690
7 1.1590000 0.2950740
8 0.1811000 0.0125870
S 8
1 2940.0000000 -0.0001230
2 441.2000000 -0.0009660
3 100.5000000 -0.0048310
4 28.4300000 -0.0193140
5 9.1690000 -0.0532800
6 3.1960000 -0.1207230
7 1.1590000 -0.1334350
8 0.1811000 0.5307670
S 1 S 1
1 0.2976000 1.0000000 1 0.0589000 1.0000000
S 1
1 0.0187700 1.0000000
P 3
1 3.6190000 0.0291110
2 0.7110000 0.1693650
3 0.1951000 0.5134580
P 1 P 1
1 1.2750000 1.0000000 1 0.0601800 1.0000000
P 1
1 0.0085000 1.0000000
D 1
1 0.2380000 1.0000000
D 1
1 0.0740000 1.0000000

View File

@ -117,18 +117,10 @@ subroutine linear_response_pp(ispin,ortho_eigvec,BSE,nBas,nC,nO,nV,nR,nOO,nVV, &
end if end if
! print*, 'pp-RPA matrix'
! call matout(nOO+nVV,nOO+nVV,M(:,:))
! write(*,*)
! Diagonalize the p-h matrix ! Diagonalize the p-h matrix
if(nOO+nVV > 0) call diagonalize_general_matrix(nOO+nVV,M,Omega,Z) if(nOO+nVV > 0) call diagonalize_general_matrix(nOO+nVV,M,Omega,Z)
! write(*,*) 'pp-RPA excitation energies'
! call matout(nOO+nVV,1,Omega(:))
! write(*,*)
! Split the various quantities in p-p and h-h parts ! Split the various quantities in p-p and h-h parts
call sort_ppRPA(ortho_eigvec,nOO,nVV,Omega(:),Z(:,:),Omega1(:),X1(:,:),Y1(:,:),Omega2(:),X2(:,:),Y2(:,:)) call sort_ppRPA(ortho_eigvec,nOO,nVV,Omega(:),Z(:,:),Omega1(:),X1(:,:),Y1(:,:),Omega2(:),X2(:,:),Y2(:,:))
@ -141,13 +133,4 @@ subroutine linear_response_pp(ispin,ortho_eigvec,BSE,nBas,nC,nO,nV,nR,nOO,nVV, &
if(abs(EcRPA - EcRPA1) > 1d-6 .or. abs(EcRPA - EcRPA2) > 1d-6) & if(abs(EcRPA - EcRPA1) > 1d-6 .or. abs(EcRPA - EcRPA2) > 1d-6) &
print*,'!!! Issue in pp-RPA linear reponse calculation RPA1 != RPA2 !!!' print*,'!!! Issue in pp-RPA linear reponse calculation RPA1 != RPA2 !!!'
! write(*,*)'X1'
! call matout(nVV,nVV,X1)
! write(*,*)'Y1'
! call matout(nVV,nOO,Y1)
! write(*,*)'X2'
! call matout(nOO,nVV,X2)
! write(*,*)'Y2'
! call matout(nOO,nOO,Y2)
end subroutine linear_response_pp end subroutine linear_response_pp

View File

@ -115,14 +115,6 @@ subroutine sort_ppRPA(ortho_eigvec,nOO,nVV,Omega,Z,Omega1,X1,Y1,Omega2,X2,Y2)
end if end if
! write(*,*) 'pp-RPA positive excitation energies'
! call matout(nVV,1,Omega1(:))
! write(*,*)
! write(*,*) 'pp-RPA negative excitation energies'
! call matout(nOO,1,Omega2(:))
! write(*,*)
! Orthogonalize eigenvectors ! Orthogonalize eigenvectors
if(ortho_eigvec) then if(ortho_eigvec) then
@ -142,8 +134,6 @@ subroutine sort_ppRPA(ortho_eigvec,nOO,nVV,Omega,Z,Omega1,X1,Y1,Omega2,X2,Y2)
ab_end = ab ab_end = ab
! print*,'deg = ',deg1,ab_start,ab_end
allocate(S1(deg1,deg1),O1(deg1,deg1)) allocate(S1(deg1,deg1),O1(deg1,deg1))
S1 = matmul(transpose(Z1(:,ab_start:ab_end)),matmul(M,Z1(:,ab_start:ab_end))) S1 = matmul(transpose(Z1(:,ab_start:ab_end)),matmul(M,Z1(:,ab_start:ab_end)))
@ -158,8 +148,6 @@ subroutine sort_ppRPA(ortho_eigvec,nOO,nVV,Omega,Z,Omega1,X1,Y1,Omega2,X2,Y2)
ab_end = ab - 1 ab_end = ab - 1
! print*,'deg = ',deg1,ab_start,ab_end
allocate(S1(deg1,deg1),O1(deg1,deg1)) allocate(S1(deg1,deg1),O1(deg1,deg1))
S1 = matmul(transpose(Z1(:,ab_start:ab_end)),matmul(M,Z1(:,ab_start:ab_end))) S1 = matmul(transpose(Z1(:,ab_start:ab_end)),matmul(M,Z1(:,ab_start:ab_end)))
@ -187,8 +175,6 @@ subroutine sort_ppRPA(ortho_eigvec,nOO,nVV,Omega,Z,Omega1,X1,Y1,Omega2,X2,Y2)
ij_end = ij ij_end = ij
! print*,'deg = ',deg2,ij_start,ij_end
allocate(S2(deg2,deg2),O2(deg2,deg2)) allocate(S2(deg2,deg2),O2(deg2,deg2))
S2 = - matmul(transpose(Z2(:,ij_start:ij_end)),matmul(M,Z2(:,ij_start:ij_end))) S2 = - matmul(transpose(Z2(:,ij_start:ij_end)),matmul(M,Z2(:,ij_start:ij_end)))
@ -203,8 +189,6 @@ subroutine sort_ppRPA(ortho_eigvec,nOO,nVV,Omega,Z,Omega1,X1,Y1,Omega2,X2,Y2)
ij_end = ij - 1 ij_end = ij - 1
! print*,'deg = ',deg2,ij_start,ij_end
allocate(S2(deg2,deg2),O2(deg2,deg2)) allocate(S2(deg2,deg2),O2(deg2,deg2))
S2 = - matmul(transpose(Z2(:,ij_start:ij_end)),matmul(M,Z2(:,ij_start:ij_end))) S2 = - matmul(transpose(Z2(:,ij_start:ij_end)),matmul(M,Z2(:,ij_start:ij_end)))
@ -239,17 +223,4 @@ subroutine sort_ppRPA(ortho_eigvec,nOO,nVV,Omega,Z,Omega1,X1,Y1,Omega2,X2,Y2)
X2(1:nVV,1:nOO) = + Z2( 1: nVV,1:nOO) X2(1:nVV,1:nOO) = + Z2( 1: nVV,1:nOO)
Y2(1:nOO,1:nOO) = - Z2(nVV+1:nOO+nVV,1:nOO) Y2(1:nOO,1:nOO) = - Z2(nVV+1:nOO+nVV,1:nOO)
! write(*,*) 'Z1t.M.Z1'
! call matout(nVV,nVV,matmul(matmul(transpose(Z1),M),Z1))
! write(*,*) 'Z2t.M.Z2'
! call matout(nOO,nOO,matmul(matmul(transpose(Z2),M),Z2))
! write(*,*) 'X1t.X1 - Y1t.Y1'
! call matout(nVV,nVV,matmul(transpose(X1),X1) - matmul(transpose(Y1),Y1))
! write(*,*) 'X2t.X2 - Y2t.Y2'
! call matout(nOO,nOO,matmul(transpose(X2),X2) - matmul(transpose(Y2),Y2))
! write(*,*) 'X1t.X2 - Y1t.Y2'
! call matout(nVV,nOO,matmul(transpose(X1),X2) - matmul(transpose(Y1),Y2))
end subroutine sort_ppRPA end subroutine sort_ppRPA