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

Accelerated selection using dichotomy

Browse Source
This commit is contained in:
Anthony Scemama 2014-06-02 15:18:45 +02:00
parent cff9652d79
commit a451a89e5a
10 changed files with 211 additions and 37 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

3
scripts/generate_h_apply.py
View File

@ -141,7 +141,7 @@ class H_apply(object):
"""
self.data["size_max"] = "256"
self.data["initialization"] = """
PROVIDE CI_electronic_energy psi_selectors_coef psi_selectors E_corr_per_selectors
PROVIDE CI_electronic_energy psi_selectors_coef psi_selectors E_corr_per_selectors psi_det_sorted_bit
"""
self.data["keys_work"] = """
call perturb_buffer_%s(i_generator,keys_out,key_idx,e_2_pert_buffer,coef_pert_buffer,sum_e_2_pert, &
@ -206,7 +206,6 @@ class H_apply(object):
call copy_h_apply_buffer_to_wf
selection_criterion_min = selection_criterion_min*0.1d0
selection_criterion = selection_criterion_min
!call remove_small_contributions
"""
self.data["keys_work"] = """
e_2_pert_buffer = 0.d0

2
src/Dets/H_apply_template.f
View File

@ -341,7 +341,7 @@ subroutine $subroutine($params_main)
$decls_main
PROVIDE H_apply_buffer_allocated mo_bielec_integrals_in_map N_det_reference psi_generators
PROVIDE H_apply_buffer_allocated mo_bielec_integrals_in_map N_det_selectors psi_generators
integer :: i_generator, k
double precision :: wall_0, wall_1, wall_2

99
src/Dets/connected_to_ref.irp.f
View File

@ -1,10 +1,77 @@
integer function connected_to_ref(key,keys,Nint,N_past_in,Ndet,thresh)
logical function is_in_wavefunction(key,keys,Nint,N_past_in,Ndet)
implicit none
integer, intent(in) :: Nint, N_past_in, Ndet
integer(bit_kind), intent(in) :: keys(Nint,2,Ndet)
integer(bit_kind), intent(in) :: key(Nint,2)
integer :: i, ibegin, iend, istep, l
integer*8 :: det_ref, det_search
integer*8, external :: det_search_key
is_in_wavefunction = .False.
ibegin = 1
iend = N_det
ASSERT (N_past > 0)
ASSERT (N_det >= N_past)
det_ref = det_search_key(key,Nint)
istep = ishft(iend-ibegin+1,-1)
i=ibegin+istep
do while (istep > 1)
i = ibegin + istep
det_search = det_search_key(psi_det_sorted_bit(1,1,i),Nint)
! print *, istep, det_ref, det_search
if ( det_search > det_ref ) then
iend = i
else if ( det_search == det_ref ) then
exit
else
ibegin = i
endif
istep = ishft(iend-ibegin+1,-1)
i = ibegin + istep
end do
! pause
do while (det_search_key(psi_det_sorted_bit(1,1,i),Nint) == det_ref)
i = i-1
if (i == 0) then
exit
endif
enddo
i += 1
do while (det_search_key(psi_det_sorted_bit(1,1,i),Nint) == det_ref)
if ( (key(1,1) /= psi_det_sorted_bit(1,1,i)).or. &
(key(1,2) /= psi_det_sorted_bit(1,2,i)) ) then
continue
else
is_in_wavefunction = .True.
!DEC$ LOOP COUNT MIN(3)
do l=2,Nint
if ( (key(l,1) /= psi_det_sorted_bit(l,1,i)).or. &
(key(l,2) /= psi_det_sorted_bit(l,2,i)) ) then
is_in_wavefunction = .False.
exit
endif
enddo
if (is_in_wavefunction) then
return
endif
endif
i += 1
enddo
end
integer function connected_to_ref(key,keys,Nint,N_past_in,Ndet)
use bitmasks
implicit none
integer, intent(in) :: Nint, N_past_in, Ndet
integer(bit_kind), intent(in) :: keys(Nint,2,Ndet)
integer(bit_kind), intent(in) :: key(Nint,2)
double precision, intent(in) :: thresh
integer :: N_past
integer :: i, l
@ -44,15 +111,21 @@ integer function connected_to_ref(key,keys,Nint,N_past_in,Ndet,thresh)
return
endif
do i=N_past,Ndet
if ( (key(1,1) /= keys(1,1,i)).or. &
(key(1,2) /= keys(1,2,i)) ) then
cycle
endif
connected_to_ref = -i
return
enddo
logical, external :: is_in_wavefunction
if (is_in_wavefunction(key,keys,Nint,N_past_in,Ndet)) then
connected_to_ref = -1
endif
return
! do i=N_past,Ndet
! if ( (key(1,1) /= keys(1,1,i)).or. &
! (key(1,2) /= keys(1,2,i)) ) then
! cycle
! endif
! connected_to_ref = -i
! return
! enddo
! return
else if (Nint==2) then
@ -80,7 +153,7 @@ integer function connected_to_ref(key,keys,Nint,N_past_in,Ndet,thresh)
endif
!DIR$ LOOP COUNT (1000)
do i=N_past+1,Ndet
do i=N_past,Ndet
if ( (key(1,1) /= keys(1,1,i)).or. &
(key(1,2) /= keys(1,2,i)).or. &
(key(2,1) /= keys(2,1,i)).or. &
@ -119,7 +192,7 @@ integer function connected_to_ref(key,keys,Nint,N_past_in,Ndet,thresh)
return
endif
do i=N_past+1,Ndet
do i=N_past,Ndet
if ( (key(1,1) /= keys(1,1,i)).or. &
(key(1,2) /= keys(1,2,i)).or. &
(key(2,1) /= keys(2,1,i)).or. &
@ -161,7 +234,7 @@ integer function connected_to_ref(key,keys,Nint,N_past_in,Ndet,thresh)
return
endif
do i=N_past+1,Ndet
do i=N_past,Ndet
if ( (key(1,1) /= keys(1,1,i)).or. &
(key(1,2) /= keys(1,2,i)) ) then
cycle

2
src/Dets/determinants.ezfio_config
View File

@ -7,4 +7,6 @@ determinants
psi_det integer*8 (determinants_N_int*determinants_bit_kind/8,2,determinants_n_det)
H_apply_threshold double precision
n_det_max_jacobi integer
threshold_generators double precision
threshold_selectors double precision

62
src/Dets/determinants.irp.f
View File

@ -90,15 +90,6 @@ END_PROVIDER
END_PROVIDER
BEGIN_PROVIDER [ integer, N_det_reference ]
implicit none
BEGIN_DOC
! Number of determinants in the reference wave function
END_DOC
N_det_reference = N_det
ASSERT (N_det_reference > 0)
END_PROVIDER
BEGIN_PROVIDER [ double precision, psi_average_norm_contrib, (N_det) ]
implicit none
BEGIN_DOC
@ -123,7 +114,7 @@ END_PROVIDER
&BEGIN_PROVIDER [ double precision, psi_average_norm_contrib_sorted, (N_det) ]
implicit none
BEGIN_DOC
! Wave function sorted by determinants (state-averaged)
! Wave function sorted by determinants contribution to the norm (state-averaged)
END_DOC
integer :: i,j,k
integer, allocatable :: iorder(:)
@ -149,3 +140,54 @@ END_PROVIDER
END_PROVIDER
BEGIN_PROVIDER [ integer(bit_kind), psi_det_sorted_bit, (N_int,2,N_det) ]
&BEGIN_PROVIDER [ double precision, psi_coef_sorted_bit, (N_det,N_states) ]
implicit none
BEGIN_DOC
! Determinants on which we apply <i|H|psi> for perturbation.
!o They are sorted by determinants interpreted as integers. Useful
! to accelerate the search of a determinant
END_DOC
integer :: i,j,k
integer, allocatable :: iorder(:)
integer*8, allocatable :: bit_tmp(:)
integer*8, external :: det_search_key
allocate ( iorder(N_det), bit_tmp(N_det) )
do i=1,N_det
iorder(i) = i
!$DIR FORCEINLINE
bit_tmp(i) = det_search_key(psi_det(1,1,i),N_int)
enddo
call i8sort(bit_tmp,iorder,N_det)
!DIR$ IVDEP
do i=1,N_det
do j=1,N_int
psi_det_sorted_bit(j,1,i) = psi_det(j,1,iorder(i))
psi_det_sorted_bit(j,2,i) = psi_det(j,2,iorder(i))
enddo
do k=1,N_states
psi_coef_sorted_bit(i,k) = psi_coef(iorder(i),k)
enddo
enddo
deallocate(iorder)
END_PROVIDER
integer*8 function det_search_key(det,Nint)
use bitmasks
implicit none
BEGIN_DOC
! Return an integer*8 corresponding to a determinant index for searching
END_DOC
integer, intent(in) :: Nint
integer(bit_kind), intent(in) :: det(Nint,2)
integer :: i
det_search_key = iand(det(1,1),det(1,2))
do i=2,Nint
det_search_key = ieor(det_search_key,iand(det(i,1),det(i,2)))
enddo
end

9
src/Full_CI/full_ci.irp.f
View File

@ -24,5 +24,14 @@ program cisd
exit
endif
enddo
print *, 'Final step'
call remove_small_contributions
call diagonalize_CI
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
deallocate(pt2,norm_pert)
end

17
src/Generators_full/generators.irp.f
View File

@ -1,14 +1,25 @@
use bitmasks
BEGIN_PROVIDER [ double precision, generators_threshold ]
BEGIN_PROVIDER [ double precision, threshold_generators ]
implicit none
BEGIN_DOC
! Percentage of the norm of the state-averaged wave function to
! consider for the generators
END_DOC
generators_threshold = 0.97d0
logical :: exists
PROVIDE ezfio_filename
call ezfio_has_determinants_threshold_generators(exists)
if (exists) then
call ezfio_get_determinants_threshold_generators(threshold_generators)
else
threshold_generators = 0.99d0
call ezfio_set_determinants_threshold_generators(threshold_generators)
endif
ASSERT (N_det > 0)
call write_double(output_Dets,threshold_generators,'Threshold on generators')
END_PROVIDER
BEGIN_PROVIDER [ integer, N_det_generators ]
implicit none
BEGIN_DOC
@ -22,7 +33,7 @@ BEGIN_PROVIDER [ integer, N_det_generators ]
N_det_generators = N_det
do i=1,N_det
norm = norm + psi_average_norm_contrib_sorted(i)
if (norm > generators_threshold) then
if (norm > threshold_generators) then
N_det_generators = i-1
exit
endif

3
src/NEEDED_MODULES
View File

@ -1 +1,2 @@
AOs Bitmask Electrons Ezfio_files MOs Nuclei Output Utils Hartree_Fock BiInts MonoInts MOGuess Dets DensityMatrix CISD Perturbation SingleRefMethod CISD_selected Selectors_full MP2 Generators_full Full_CI CISD_SC2_selected CISD_SC2_selected
AOs BiInts Bitmask CISD CISD_SC2_selected CISD_selected DensityMatrix Dets Electrons Ezfio_files Full_CI Generators_full Hartree_Fock MOGuess MonoInts MOs MP2 Nuclei Output Perturbation SC2 Selectors_full SingleRefMethod Utils

4
src/Perturbation/perturbation_template.f
View File

@ -15,7 +15,7 @@ subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,c
double precision, intent(inout) :: coef_pert_buffer(N_st,buffer_size),e_2_pert_buffer(N_st,buffer_size),sum_H_pert_diag(N_st)
double precision :: c_pert(N_st), e_2_pert(N_st), H_pert_diag(N_st)
integer :: i,k, c_ref
integer :: connected_to_ref
integer, external :: connected_to_ref
ASSERT (Nint > 0)
ASSERT (Nint == N_int)
@ -24,7 +24,7 @@ subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,c
ASSERT (N_st > 0)
do i = 1,buffer_size
c_ref = connected_to_ref(buffer(1,1,i),psi_generators,Nint,i_generator,N_det,h_apply_threshold)
c_ref = connected_to_ref(buffer(1,1,i),psi_generators,Nint,i_generator,N_det)
if (c_ref /= 0) then
cycle

47
src/Selectors_full/selectors.irp.f
View File

@ -1,31 +1,68 @@
use bitmasks
BEGIN_PROVIDER [ double precision, threshold_selectors ]
implicit none
BEGIN_DOC
! Percentage of the norm of the state-averaged wave function to
! consider for the selectors
END_DOC
logical :: exists
PROVIDE ezfio_filename
call ezfio_has_determinants_threshold_selectors(exists)
if (exists) then
call ezfio_get_determinants_threshold_selectors(threshold_selectors)
else
threshold_selectors = 0.99d0
call ezfio_set_determinants_threshold_selectors(threshold_selectors)
endif
ASSERT (N_det > 0)
call write_double(output_Dets,threshold_selectors,'Threshold on selectors')
END_PROVIDER
BEGIN_PROVIDER [ integer, psi_selectors_size ]
implicit none
psi_selectors_size = psi_det_size
END_PROVIDER
BEGIN_PROVIDER [ integer, N_det_selectors]
implicit none
N_det_selectors = N_det
BEGIN_DOC
! For Single reference wave functions, the number of selectors is 1 : the
! Hartree-Fock determinant
END_DOC
integer :: i
double precision :: norm
call write_time(output_dets)
norm = 0.d0
N_det_selectors = N_det
do i=1,N_det
norm = norm + psi_average_norm_contrib_sorted(i)
if (norm > threshold_selectors) then
N_det_selectors = i-1
exit
endif
enddo
N_det_selectors = max(N_det_selectors,1)
call write_int(output_dets,N_det_selectors,'Number of selectors')
END_PROVIDER
BEGIN_PROVIDER [ integer(bit_kind), psi_selectors, (N_int,2,psi_selectors_size) ]
&BEGIN_PROVIDER [ double precision, psi_selectors_coef, (psi_selectors_size,N_states) ]
implicit none
BEGIN_DOC
! On what we apply <i|H|psi> for perturbation. If selection, it may be 0.9 of the norm.
! Determinants on which we apply <i|H|psi> for perturbation.
END_DOC
integer :: i,k
do i=1,N_det_selectors
do k=1,N_int
psi_selectors(k,1,i) = psi_det(k,1,i)
psi_selectors(k,2,i) = psi_det(k,2,i)
psi_selectors(k,1,i) = psi_det_sorted(k,1,i)
psi_selectors(k,2,i) = psi_det_sorted(k,2,i)
enddo
enddo
do k=1,N_states
do i=1,N_det_selectors
psi_selectors_coef(i,k) = psi_coef(i,k)
psi_selectors_coef(i,k) = psi_coef_sorted(i,k)
enddo
enddo
END_PROVIDER