10
0
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-08-24 21:21:47 +02:00

added print_hmat and sparse_mat

This commit is contained in:
eginer 2022-12-21 13:54:20 +01:00
parent a3bc5fd421
commit f7e58e4a63
5 changed files with 279 additions and 1 deletions

View File

@ -96,7 +96,6 @@ subroutine filter_not_connected(key1,key2,Nint,sze,idx)
idx(0) = l-1
end
subroutine filter_connected(key1,key2,Nint,sze,idx)
use bitmasks
implicit none

View File

@ -0,0 +1,164 @@
use bitmasks
subroutine filter_connected_array(key1,key2,ld,Nint,sze,idx)
use bitmasks
implicit none
BEGIN_DOC
! Filters out the determinants that are not connected by H
!
! returns the array idx which contains the index of the
!
! determinants in the array key1 that interact
!
! via the H operator with key2.
!
! idx(0) is the number of determinants that interact with key1
END_DOC
integer, intent(in) :: Nint, ld,sze
integer(bit_kind), intent(in) :: key1(Nint,2,ld)
integer(bit_kind), intent(in) :: key2(Nint,2)
integer, intent(out) :: idx(0:sze)
integer :: i,j,l
integer :: degree_x2
ASSERT (Nint > 0)
ASSERT (sze >= 0)
l=1
if (Nint==1) then
!DIR$ LOOP COUNT (1000)
do i=1,sze
degree_x2 = popcnt( xor( key1(1,1,i), key2(1,1))) &
+ popcnt( xor( key1(1,2,i), key2(1,2)))
! print*,degree_x2
if (degree_x2 > 4) then
cycle
else
idx(l) = i
l = l+1
endif
enddo
else if (Nint==2) then
!DIR$ LOOP COUNT (1000)
do i=1,sze
degree_x2 = popcnt(xor( key1(1,1,i), key2(1,1))) + &
popcnt(xor( key1(2,1,i), key2(2,1))) + &
popcnt(xor( key1(1,2,i), key2(1,2))) + &
popcnt(xor( key1(2,2,i), key2(2,2)))
if (degree_x2 > 4) then
cycle
else
idx(l) = i
l = l+1
endif
enddo
else if (Nint==3) then
!DIR$ LOOP COUNT (1000)
do i=1,sze
degree_x2 = 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)))
if (degree_x2 > 4) then
cycle
else
idx(l) = i
l = l+1
endif
enddo
else
!DIR$ LOOP COUNT (1000)
do i=1,sze
degree_x2 = 0
!DIR$ LOOP COUNT MIN(4)
do j=1,Nint
degree_x2 = degree_x2+ popcnt(xor( key1(j,1,i), key2(j,1))) +&
popcnt(xor( key1(j,2,i), key2(j,2)))
if (degree_x2 > 4) then
exit
endif
enddo
if (degree_x2 <= 5) then
idx(l) = i
l = l+1
endif
enddo
endif
idx(0) = l-1
! print*,'idx(0) = ',idx(0)
end
BEGIN_PROVIDER [ integer, n_sparse_mat]
&BEGIN_PROVIDER [ integer, n_connected_per_det, (N_det)]
&BEGIN_PROVIDER [ integer, n_max_connected_per_det]
implicit none
BEGIN_DOC
! n_sparse_mat = total number of connections in the CI matrix
!
! n_connected_per_det(i) = number of connected determinants to the determinant psi_det(1,1,i)
!
! n_max_connected_per_det = maximum number of connected determinants
END_DOC
integer, allocatable :: idx(:)
allocate(idx(0:N_det))
integer :: i
n_sparse_mat = 0
do i = 1, N_det
call filter_connected_array(psi_det_sorted,psi_det_sorted(1,1,i),psi_det_size,N_int,N_det,idx)
n_connected_per_det(i) = idx(0)
n_sparse_mat += idx(0)
enddo
n_max_connected_per_det = maxval(n_connected_per_det)
END_PROVIDER
BEGIN_PROVIDER [ integer(bit_kind), connected_det_per_det, (N_int,2,n_max_connected_per_det,N_det)]
&BEGIN_PROVIDER [ integer(bit_kind), list_connected_det_per_det, (n_max_connected_per_det,N_det)]
implicit none
BEGIN_DOC
! connected_det_per_det(:,:,j,i) = jth connected determinant to the determinant psi_det(:,:,i)
!
! list_connected_det_per_det(j,i) = index of jth determinant in psi_det which is connected to psi_det(:,:,i)
END_DOC
integer, allocatable :: idx(:)
allocate(idx(0:N_det))
integer :: i,j
do i = 1, N_det
call filter_connected_array(psi_det_sorted,psi_det_sorted(1,1,i),psi_det_size,N_int,N_det,idx)
do j = 1, idx(0)
connected_det_per_det(1:N_int,1:2,j,i) = psi_det_sorted(1:N_int,1:2,idx(j))
list_connected_det_per_det(j,i) = idx(j)
enddo
enddo
END_PROVIDER
BEGIN_PROVIDER [ double precision, sparse_h_mat, (n_max_connected_per_det, N_det)]
implicit none
BEGIN_DOC
! sparse matrix format
!
! sparse_h_mat(j,i) = matrix element between the jth connected determinant and psi_det(:,:,i)
END_DOC
integer :: i,j
double precision :: hij
do i = 1, N_det
do j = 1, n_connected_per_det(i)
call i_H_j(psi_det(1,1,i),connected_det_per_det(1,1,j,i),N_int,hij)
sparse_h_mat(j,i) = hij
enddo
enddo
END_PROVIDER

View File

@ -73,6 +73,29 @@
+ (Fock_matrix_tc_mo_beta(i,j) - Fock_matrix_tc_mo_alpha(i,j))
enddo
enddo
if(three_body_h_tc)then
! C-O
do j = 1, elec_beta_num
do i = elec_beta_num+1, elec_alpha_num
Fock_matrix_tc_mo_tot(i,j) += 0.5d0*(fock_a_tot_3e_bi_orth(i,j) + fock_b_tot_3e_bi_orth(i,j))
Fock_matrix_tc_mo_tot(j,i) += 0.5d0*(fock_a_tot_3e_bi_orth(j,i) + fock_b_tot_3e_bi_orth(j,i))
enddo
enddo
! C-V
do j = 1, elec_beta_num
do i = elec_alpha_num+1, mo_num
Fock_matrix_tc_mo_tot(i,j) += 0.5d0*(fock_a_tot_3e_bi_orth(i,j) + fock_b_tot_3e_bi_orth(i,j))
Fock_matrix_tc_mo_tot(j,i) += 0.5d0*(fock_a_tot_3e_bi_orth(j,i) + fock_b_tot_3e_bi_orth(j,i))
enddo
enddo
! O-V
do j = elec_beta_num+1, elec_alpha_num
do i = elec_alpha_num+1, mo_num
Fock_matrix_tc_mo_tot(i,j) += 0.5d0*(fock_a_tot_3e_bi_orth(i,j) + fock_b_tot_3e_bi_orth(i,j))
Fock_matrix_tc_mo_tot(j,i) += 0.5d0*(fock_a_tot_3e_bi_orth(j,i) + fock_b_tot_3e_bi_orth(j,i))
enddo
enddo
endif
endif

View File

@ -128,6 +128,8 @@ BEGIN_PROVIDER [double precision, diag_three_elem_hf]
call give_abb_contrib(integral_abb)
call give_bbb_contrib(integral_bbb)
diag_three_elem_hf = integral_aaa + integral_aab + integral_abb + integral_bbb
! print*,'integral_aaa + integral_aab + integral_abb + integral_bbb'
! print*,integral_aaa , integral_aab , integral_abb , integral_bbb
endif

View File

@ -0,0 +1,90 @@
program print_h_mat
implicit none
BEGIN_DOC
! program that prints out the CI matrix in sparse form
END_DOC
read_wf = .True.
touch read_wf
call print_wf_dets
call print_wf_coef
call sparse_mat
call full_mat
call test_sparse_mat
end
subroutine print_wf_dets
implicit none
integer :: i,j
character*(128) :: output
integer :: i_unit_output,getUnitAndOpen
output=trim(ezfio_filename)//'.wf_det'
i_unit_output = getUnitAndOpen(output,'w')
write(i_unit_output,*)N_det,N_int
do i = 1, N_det
write(i_unit_output,*)psi_det_sorted(1:N_int,1,i)
write(i_unit_output,*)psi_det_sorted(1:N_int,2,i)
enddo
end
subroutine print_wf_coef
implicit none
integer :: i,j
character*(128) :: output
integer :: i_unit_output,getUnitAndOpen
output=trim(ezfio_filename)//'.wf_coef'
i_unit_output = getUnitAndOpen(output,'w')
write(i_unit_output,*)N_det,N_states
do i = 1, N_det
write(i_unit_output,*)psi_coef_sorted(i,1:N_states)
enddo
end
subroutine sparse_mat
implicit none
integer :: i,j
character*(128) :: output
integer :: i_unit_output,getUnitAndOpen
output=trim(ezfio_filename)//'.hmat_sparse'
i_unit_output = getUnitAndOpen(output,'w')
do i = 1, N_det
write(i_unit_output,*)i,n_connected_per_det(i)
do j =1, n_connected_per_det(i)
write(i_unit_output,*)list_connected_det_per_det(j,i),sparse_h_mat(j,i)
enddo
enddo
end
subroutine full_mat
implicit none
integer :: i,j
character*(128) :: output
integer :: i_unit_output,getUnitAndOpen
output=trim(ezfio_filename)//'.hmat_full'
i_unit_output = getUnitAndOpen(output,'w')
do i = 1, N_det
do j = i, N_det
write(i_unit_output,*)i,j,H_matrix_all_dets(j,i)
enddo
enddo
end
subroutine test_sparse_mat
implicit none
integer :: i,j
double precision, allocatable :: eigvec(:,:), eigval(:), hmat(:,:)
allocate(eigval(N_det), eigvec(N_det,N_det),hmat(N_det,N_det))
hmat = 0.d0
do i = 1, N_det
do j =1, n_connected_per_det(i)
hmat(list_connected_det_per_det(j,i),i) = sparse_h_mat(j,i)
enddo
enddo
call lapack_diag(eigval,eigvec,hmat,N_det,N_det)
print*,'The two energies should be the same '
print*,'eigval(1) = ',eigval(1)
print*,'psi_energy= ',CI_electronic_energy(1)
end