From 99aa50048396e09204739a1c9c1d84068abfc0bd Mon Sep 17 00:00:00 2001 From: Pierre-Francois Loos Date: Mon, 28 Nov 2022 15:21:05 +0100 Subject: [PATCH] missing files for BSE2@GW --- src/GW/BSE2_static_kernel_KA.f90 | 91 ++++++++++++++++++++++++++++++++ src/GW/BSE2_static_kernel_KB.f90 | 91 ++++++++++++++++++++++++++++++++ src/GW/static_screening_W.f90 | 58 ++++++++++++++++++++ 3 files changed, 240 insertions(+) create mode 100644 src/GW/BSE2_static_kernel_KA.f90 create mode 100644 src/GW/BSE2_static_kernel_KB.f90 create mode 100644 src/GW/static_screening_W.f90 diff --git a/src/GW/BSE2_static_kernel_KA.f90 b/src/GW/BSE2_static_kernel_KA.f90 new file mode 100644 index 0000000..2f8fec5 --- /dev/null +++ b/src/GW/BSE2_static_kernel_KA.f90 @@ -0,0 +1,91 @@ +subroutine BSE2_static_kernel_KA(eta,nBas,nC,nO,nV,nR,nS,lambda,eW,W,KA2_sta) + +! Compute the second-order static BSE kernel for the resonant block (only for singlets!) + + implicit none + include 'parameters.h' + +! Input variables + + double precision,intent(in) :: eta + integer,intent(in) :: nBas + integer,intent(in) :: nC + integer,intent(in) :: nO + integer,intent(in) :: nV + integer,intent(in) :: nR + integer,intent(in) :: nS + double precision,intent(in) :: lambda + + double precision,intent(in) :: eW(nBas) + double precision,intent(in) :: W(nBas,nBas,nBas,nBas) + +! Local variables + + double precision :: dem,num + integer :: i,j,k,l + integer :: a,b,c,d + integer :: ia,jb + +! Output variables + + double precision,intent(out) :: KA2_sta(nS,nS) + +!------------------------------------------------ +! Compute BSE2 kernel +!------------------------------------------------ + + ia = 0 + do i=nC+1,nO + do a=nO+1,nBas-nR + ia = ia + 1 + + jb = 0 + do j=nC+1,nO + do b=nO+1,nBas-nR + jb = jb + 1 + + do k=nC+1,nO + do c=nO+1,nBas-nR + + dem = - (eW(c) - eW(k)) + num = 2d0*W(j,k,i,c)*W(a,c,b,k) + + KA2_sta(ia,jb) = KA2_sta(ia,jb) - num*dem/(dem**2 + eta**2) + + dem = + (eW(c) - eW(k)) + num = 2d0*W(j,c,i,k)*W(a,k,b,c) + + KA2_sta(ia,jb) = KA2_sta(ia,jb) + num*dem/(dem**2 + eta**2) + + end do + end do + + do c=nO+1,nBas-nR + do d=nO+1,nBas-nR + + dem = - (eW(c) + eW(d)) + num = 2d0*W(a,j,c,d)*W(c,d,i,b) + + KA2_sta(ia,jb) = KA2_sta(ia,jb) + num*dem/(dem**2 + eta**2) + + end do + end do + + do k=nC+1,nO + do l=nC+1,nO + + dem = - (eW(k) + eW(l)) + num = 2d0*W(a,j,k,l)*W(k,l,i,b) + + KA2_sta(ia,jb) = KA2_sta(ia,jb) - num*dem/(dem**2 + eta**2) + + end do + end do + + end do + end do + + end do + end do + +end subroutine BSE2_static_kernel_KA diff --git a/src/GW/BSE2_static_kernel_KB.f90 b/src/GW/BSE2_static_kernel_KB.f90 new file mode 100644 index 0000000..cb4abe3 --- /dev/null +++ b/src/GW/BSE2_static_kernel_KB.f90 @@ -0,0 +1,91 @@ +subroutine BSE2_static_kernel_KB(eta,nBas,nC,nO,nV,nR,nS,lambda,eW,W,KB2_sta) + +! Compute the second-order static BSE kernel for the coupling block (only for singlets!) + + implicit none + include 'parameters.h' + +! Input variables + + double precision,intent(in) :: eta + integer,intent(in) :: nBas + integer,intent(in) :: nC + integer,intent(in) :: nO + integer,intent(in) :: nV + integer,intent(in) :: nR + integer,intent(in) :: nS + double precision,intent(in) :: lambda + + double precision,intent(in) :: eW(nBas) + double precision,intent(in) :: W(nBas,nBas,nBas,nBas) + +! Local variables + + double precision :: dem,num + integer :: i,j,k,l + integer :: a,b,c,d + integer :: ia,jb + +! Output variables + + double precision,intent(out) :: KB2_sta(nS,nS) + +!------------------------------------------------ +! Compute BSE2 kernel +!------------------------------------------------ + + ia = 0 + do i=nC+1,nO + do a=nO+1,nBas-nR + ia = ia + 1 + + jb = 0 + do j=nC+1,nO + do b=nO+1,nBas-nR + jb = jb + 1 + + do k=nC+1,nO + do c=nO+1,nBas-nR + + dem = - (eW(c) - eW(k)) + num = 2d0*W(b,k,i,c)*W(a,c,j,k) + + KB2_sta(ia,jb) = KB2_sta(ia,jb) - num*dem/(dem**2 + eta**2) + + dem = + (eW(c) - eW(k)) + num = 2d0*W(b,c,i,k)*W(a,k,b,c) + + KB2_sta(ia,jb) = KB2_sta(ia,jb) + num*dem/(dem**2 + eta**2) + + end do + end do + + do c=nO+1,nBas-nR + do d=nO+1,nBas-nR + + dem = - (eW(c) + eW(d)) + num = 2d0*W(a,b,c,d)*W(c,d,i,j) + + KB2_sta(ia,jb) = KB2_sta(ia,jb) + num*dem/(dem**2 + eta**2) + + end do + end do + + do k=nC+1,nO + do l=nC+1,nO + + dem = - (eW(k) + eW(l)) + num = 2d0*W(a,b,k,l)*W(k,l,i,j) + + KB2_sta(ia,jb) = KB2_sta(ia,jb) - num*dem/(dem**2 + eta**2) + + end do + end do + + end do + end do + + end do + end do + +end subroutine BSE2_static_kernel_KB diff --git a/src/GW/static_screening_W.f90 b/src/GW/static_screening_W.f90 new file mode 100644 index 0000000..9d9942e --- /dev/null +++ b/src/GW/static_screening_W.f90 @@ -0,0 +1,58 @@ +subroutine static_kernel_W(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,Om,rho,W) + +! Compute the second-order static BSE kernel for the resonant block (only for singlets!) + + implicit none + include 'parameters.h' + +! Input variables + + double precision,intent(in) :: eta + integer,intent(in) :: nBas + integer,intent(in) :: nC + integer,intent(in) :: nO + integer,intent(in) :: nV + integer,intent(in) :: nR + integer,intent(in) :: nS + double precision,intent(in) :: lambda + + double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas) + double precision,intent(in) :: Om(nS) + + double precision,intent(in) :: rho(nBas,nBas,nS) + + +! Local variables + + double precision :: chi + integer :: p,q,r,s + integer :: m + double precision :: dem + +! Output variables + + double precision,intent(out) :: W(nBas,nBas,nBas,nBas) + +!------------------------------------------------ +! Compute static screening (physicist's notation) +!------------------------------------------------ + + do p=1,nBas + do q=1,nBas + do r=1,nBas + do s=1,nBas + + chi = 0d0 + do m=1,nS + dem = Om(m)**2 + eta**2 + chi = chi + rho(p,q,m)*rho(r,s,m)*Om(m)/dem + enddo + + W(p,s,q,r) = - lambda*ERI(p,s,q,r) + 4d0*lambda*chi + + enddo + enddo + enddo + enddo + +end subroutine static_kernel_W