LCPQ/quantum_package
10
0
Fork 0
mirror of https://github.com/LCPQ/quantum_package synced 2024-06-22 05:02:15 +02:00
Code Issues Releases Wiki Activity

Accelerated get_excitation_degree_vector

Browse Source
This commit is contained in:
Anthony Scemama 2015-07-29 15:35:01 +02:00
parent c042e5af0c
commit 7ef8ee171d
2 changed files with 118 additions and 98 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

170
plugins/MRCC_Utils_new/mrcc_amplitudes.irp.f
View File

@ -1,85 +1,93 @@
subroutine get_excitation_operators_for_one_ref(det_ref,i_state,ndetnonref,N_connect_ref,excitation_operators,amplitudes_phase_less,index_connected)
use bitmasks
implicit none
integer(bit_kind), intent(in) :: det_ref(N_int,2)
integer, intent(in) :: i_state,ndetnonref
integer*2, intent(out) :: excitation_operators(5,ndetnonref)
integer, intent(out) :: index_connected(ndetnonref)
integer, intent(out) :: N_connect_ref
double precision, intent(out) :: amplitudes_phase_less(ndetnonref)
integer :: i,j,k,l,degree,h1,p1,h2,p2,s1,s2
integer :: exc(0:2,2,2)
double precision :: phase,hij
BEGIN_DOC
! This subroutine provides all the amplitudes and excitation operators
! that one needs to go from the reference to the non reference wave function
! you enter with det_ref that is a reference determinant
!
! N_connect_ref is the number of determinants belonging to psi_non_ref
! that are connected to det_ref.
!
! amplitudes_phase_less(i) = amplitude phase less t_{I->i} = <I|H|i> * lambda_mrcc(i) * phase(I->i)
!
! excitation_operators(:,i) represents the holes and particles that
! link the ith connected determinant to det_ref
! if ::
! excitation_operators(5,i) = 2 :: double excitation alpha
! excitation_operators(5,i) = -2 :: double excitation beta
!!! excitation_operators(1,i) :: hole 1
!!! excitation_operators(2,i) :: particle 1
!!! excitation_operators(3,i) :: hole 2
!!! excitation_operators(4,i) :: particle 2
! else if ::
! excitation_operators(5,i) = 1 :: single excitation alpha
!!! excitation_operators(1,i) :: hole 1
!!! excitation_operators(2,i) :: particle 1
! else if ::
! excitation_operators(5,i) = -1 :: single excitation beta
!!! excitation_operators(3,i) :: hole 1
!!! excitation_operators(4,i) :: particle 1
! else if ::
!!! excitation_operators(5,i) = 0 :: double excitation alpha/beta
!!! excitation_operators(1,i) :: hole 1 alpha
!!! excitation_operators(2,i) :: particle 1 alpha
!!! excitation_operators(3,i) :: hole 2 beta
!!! excitation_operators(4,i) :: particle 2 beta
END_DOC
N_connect_ref = 0
do i = 1, ndetnonref
call i_H_j_phase_out(det_ref,psi_non_ref(1,1,i),N_int,hij,phase,exc,degree)
! if(dabs(hij).le.mo_integrals_threshold)cycle
N_connect_ref +=1
index_connected(N_connect_ref) = i
call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2)
amplitudes_phase_less(N_connect_ref) = hij * lambda_mrcc(i_state,i) !*phase
if(degree==2)then
excitation_operators(1,N_connect_ref) = h1
excitation_operators(2,N_connect_ref) = p1
excitation_operators(3,N_connect_ref) = h2
excitation_operators(4,N_connect_ref) = p2
if(s1==s2.and.s1==1)then ! double alpha
excitation_operators(5,N_connect_ref)= 2
elseif(s1==s2.and.s1==2)then ! double beta
excitation_operators(5,N_connect_ref)=-2
else
excitation_operators(5,N_connect_ref)= 0 ! double alpha/beta
endif
elseif(degree==1)then
if(s1==1)then ! mono alpha
excitation_operators(5,N_connect_ref)= 1
excitation_operators(1,N_connect_ref) = h1
excitation_operators(2,N_connect_ref) = p1
else ! mono beta
excitation_operators(5,N_connect_ref)=-1
excitation_operators(3,N_connect_ref) = h1
excitation_operators(4,N_connect_ref) = p1
endif
else
N_connect_ref-=1
endif
use bitmasks
implicit none
integer(bit_kind), intent(in) :: det_ref(N_int,2)
integer, intent(in) :: i_state,ndetnonref
integer*2, intent(out) :: excitation_operators(5,ndetnonref)
integer, intent(out) :: index_connected(ndetnonref)
integer, intent(out) :: N_connect_ref
double precision, intent(out) :: amplitudes_phase_less(ndetnonref)
enddo
integer :: i,j,k,l,degree,h1,p1,h2,p2,s1,s2
integer :: exc(0:2,2,2)
double precision :: phase,hij
BEGIN_DOC
! This subroutine provides all the amplitudes and excitation operators
! that one needs to go from the reference to the non reference wave function
! you enter with det_ref that is a reference determinant
!
! N_connect_ref is the number of determinants belonging to psi_non_ref
! that are connected to det_ref.
!
! amplitudes_phase_less(i) = amplitude phase less t_{I->i} = <I|H|i> * lambda_mrcc(i) * phase(I->i)
!
! excitation_operators(:,i) represents the holes and particles that
! link the ith connected determinant to det_ref
! if ::
! excitation_operators(5,i) = 2 :: double excitation alpha
! excitation_operators(5,i) = -2 :: double excitation beta
!!! excitation_operators(1,i) :: hole 1
!!! excitation_operators(2,i) :: particle 1
!!! excitation_operators(3,i) :: hole 2
!!! excitation_operators(4,i) :: particle 2
! else if ::
! excitation_operators(5,i) = 1 :: single excitation alpha
!!! excitation_operators(1,i) :: hole 1
!!! excitation_operators(2,i) :: particle 1
! else if ::
! excitation_operators(5,i) = -1 :: single excitation beta
!!! excitation_operators(3,i) :: hole 1
!!! excitation_operators(4,i) :: particle 1
! else if ::
!!! excitation_operators(5,i) = 0 :: double excitation alpha/beta
!!! excitation_operators(1,i) :: hole 1 alpha
!!! excitation_operators(2,i) :: particle 1 alpha
!!! excitation_operators(3,i) :: hole 2 beta
!!! excitation_operators(4,i) :: particle 2 beta
END_DOC
N_connect_ref = 0
do i = 1, ndetnonref
call i_H_j_phase_out(det_ref,psi_non_ref(1,1,i),N_int,hij,phase,exc,degree)
if (dabs(hij) <= mo_integrals_threshold) then
cycle
endif
N_connect_ref +=1
index_connected(N_connect_ref) = i
call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2)
amplitudes_phase_less(N_connect_ref) = hij * lambda_mrcc(i_state,i) !*phase
if (degree==2) then
excitation_operators(1,N_connect_ref) = h1
excitation_operators(2,N_connect_ref) = p1
excitation_operators(3,N_connect_ref) = h2
excitation_operators(4,N_connect_ref) = p2
if(s1==s2.and.s1==1)then ! double alpha
excitation_operators(5,N_connect_ref) = 2
elseif(s1==s2.and.s1==2)then ! double beta
excitation_operators(5,N_connect_ref) = -2
else ! double alpha/beta
excitation_operators(5,N_connect_ref) = 0
endif
else if(degree==1) then
if(s1==1)then ! mono alpha
excitation_operators(5,N_connect_ref) = 1
excitation_operators(1,N_connect_ref) = h1
excitation_operators(2,N_connect_ref) = p1
else ! mono beta
excitation_operators(5,N_connect_ref) = -1
excitation_operators(3,N_connect_ref) = h1
excitation_operators(4,N_connect_ref) = p1
endif
else
N_connect_ref-=1
endif
enddo
end

46
src/Determinants/slater_rules.irp.f
View File

@ -960,7 +960,7 @@ subroutine get_excitation_degree_vector(key1,key2,degree,Nint,sze,idx)
integer, intent(out) :: degree(sze)
integer, intent(out) :: idx(0:sze)
integer :: i,l
integer :: i,l,d
ASSERT (Nint > 0)
ASSERT (sze > 0)
@ -970,9 +970,12 @@ subroutine get_excitation_degree_vector(key1,key2,degree,Nint,sze,idx)
!DIR$ LOOP COUNT (1000)
do i=1,sze
degree(l) = ishft(popcnt(xor( key1(1,1,i), key2(1,1))) + &
d = ishft(popcnt(xor( key1(1,1,i), key2(1,1))) + &
popcnt(xor( key1(1,2,i), key2(1,2))),-1)
if (degree(l) < 3) then
if (d > 2) then
cycle
else
degree(l) = d
idx(l) = i
l = l+1
endif
@ -982,13 +985,16 @@ subroutine get_excitation_degree_vector(key1,key2,degree,Nint,sze,idx)
!DIR$ LOOP COUNT (1000)
do i=1,sze
degree(l) = ishft(popcnt(xor( key1(1,1,i), key2(1,1))) + &
d = ishft(popcnt(xor( key1(1,1,i), key2(1,1))) + &
popcnt(xor( key1(1,2,i), key2(1,2))) + &
popcnt(xor( key1(2,1,i), key2(2,1))) + &
popcnt(xor( key1(2,2,i), key2(2,2))),-1)
if (degree(l) < 3) then
idx(l) = i
l = l+1
if (d > 2) then
cycle
else
degree(l) = d
idx(l) = i
l = l+1
endif
enddo
@ -996,15 +1002,18 @@ subroutine get_excitation_degree_vector(key1,key2,degree,Nint,sze,idx)
!DIR$ LOOP COUNT (1000)
do i=1,sze
degree(l) = ishft( popcnt(xor( key1(1,1,i), key2(1,1))) + &
d = ishft( popcnt(xor( key1(1,1,i), key2(1,1))) + &
popcnt(xor( key1(1,2,i), key2(1,2))) + &
popcnt(xor( key1(2,1,i), key2(2,1))) + &
popcnt(xor( key1(2,2,i), key2(2,2))) + &
popcnt(xor( key1(3,1,i), key2(3,1))) + &
popcnt(xor( key1(3,2,i), key2(3,2))),-1)
if (degree(l) < 3) then
idx(l) = i
l = l+1
if (d > 2) then
cycle
else
degree(l) = d
idx(l) = i
l = l+1
endif
enddo
@ -1012,16 +1021,19 @@ subroutine get_excitation_degree_vector(key1,key2,degree,Nint,sze,idx)
!DIR$ LOOP COUNT (1000)
do i=1,sze
degree(l) = 0
d = 0
!DEC$ LOOP COUNT MIN(4)
do l=1,Nint
degree(l) = degree(l)+ popcnt(xor( key1(l,1,i), key2(l,1))) +&
d = d + popcnt(xor( key1(l,1,i), key2(l,1))) +&
popcnt(xor( key1(l,2,i), key2(l,2)))
enddo
degree(l) = ishft(degree(l),-1)
if (degree(l) < 3) then
idx(l) = i
l = l+1
d = ishft(d,-1)
if (d > 2) then
cycle
else
degree(l) = d
idx(l) = i
l = l+1
endif
enddo