diff --git a/config/gfortran.cfg b/config/gfortran.cfg index 066639a8..a4bca336 100644 --- a/config/gfortran.cfg +++ b/config/gfortran.cfg @@ -10,7 +10,7 @@ # # [COMMON] -FC : gfortran -g -ffree-line-length-none -I . -static-libgcc +FC : gfortran -g -ffree-line-length-none -I . -static-libgcc -msse4_2 LAPACK_LIB : -llapack -lblas IRPF90 : irpf90 IRPF90_FLAGS : --ninja --align=32 @@ -22,7 +22,7 @@ IRPF90_FLAGS : --ninja --align=32 # 0 : Deactivate # [OPTION] -MODE : OPT ; [ OPT | PROFILE | DEBUG ] : Chooses the section below +MODE : DEBUG ; [ OPT | PROFILE | DEBUG ] : Chooses the section below CACHE : 1 ; Enable cache_compile.py OPENMP : 1 ; Append OpenMP flags diff --git a/plugins/MRCC_Utils/mrcc_dress.irp.f b/plugins/MRCC_Utils/mrcc_dress.irp.f index a5f9e068..42c8d9d0 100644 --- a/plugins/MRCC_Utils/mrcc_dress.irp.f +++ b/plugins/MRCC_Utils/mrcc_dress.irp.f @@ -14,52 +14,52 @@ BEGIN_PROVIDER [ integer(omp_lock_kind), psi_ref_lock, (psi_det_size) ] END_PROVIDER -subroutine create_minilist(key_mask, fullList, miniList, idx_miniList, N_fullList, N_miniList, Nint) - use bitmasks - implicit none - - integer(bit_kind), intent(in) :: fullList(Nint, 2, N_fullList) - integer, intent(in) :: N_fullList - integer(bit_kind),intent(out) :: miniList(Nint, 2, N_fullList) - integer,intent(out) :: idx_miniList(N_fullList), N_miniList - integer, intent(in) :: Nint - integer(bit_kind) :: key_mask(Nint, 2) - integer :: ni, i, n_a, n_b, e_a, e_b - - - n_a = 0 - n_b = 0 - do ni=1,nint - n_a = n_a + popcnt(key_mask(ni,1)) - n_b = n_b + popcnt(key_mask(ni,2)) - end do - - if(n_a == 0) then - N_miniList = N_fullList - miniList(:,:,:) = fullList(:,:,:) - do i=1,N_fullList - idx_miniList(i) = i - end do - return - end if - - N_miniList = 0 - - do i=1,N_fullList - e_a = n_a - e_b = n_b - do ni=1,nint - e_a -= popcnt(iand(fullList(ni, 1, i), key_mask(ni, 1))) - e_b -= popcnt(iand(fullList(ni, 2, i), key_mask(ni, 2))) - end do - - if(e_a + e_b <= 2) then - N_miniList = N_miniList + 1 - miniList(:,:,N_miniList) = fullList(:,:,i) - idx_miniList(N_miniList) = i - end if - end do -end subroutine +! subroutine create_minilist(key_mask, fullList, miniList, idx_miniList, N_fullList, N_miniList, Nint) +! use bitmasks +! implicit none +! +! integer(bit_kind), intent(in) :: fullList(Nint, 2, N_fullList) +! integer, intent(in) :: N_fullList +! integer(bit_kind),intent(out) :: miniList(Nint, 2, N_fullList) +! integer,intent(out) :: idx_miniList(N_fullList), N_miniList +! integer, intent(in) :: Nint +! integer(bit_kind) :: key_mask(Nint, 2) +! integer :: ni, i, n_a, n_b, e_a, e_b +! +! +! n_a = 0 +! n_b = 0 +! do ni=1,nint +! n_a = n_a + popcnt(key_mask(ni,1)) +! n_b = n_b + popcnt(key_mask(ni,2)) +! end do +! +! if(n_a == 0) then +! N_miniList = N_fullList +! miniList(:,:,:) = fullList(:,:,:) +! do i=1,N_fullList +! idx_miniList(i) = i +! end do +! return +! end if +! +! N_miniList = 0 +! +! do i=1,N_fullList +! e_a = n_a +! e_b = n_b +! do ni=1,nint +! e_a -= popcnt(iand(fullList(ni, 1, i), key_mask(ni, 1))) +! e_b -= popcnt(iand(fullList(ni, 2, i), key_mask(ni, 2))) +! end do +! +! if(e_a + e_b <= 2) then +! N_miniList = N_miniList + 1 +! miniList(:,:,N_miniList) = fullList(:,:,i) +! idx_miniList(N_miniList) = i +! end if +! end do +! end subroutine subroutine mrcc_dress(delta_ij_, delta_ii_, Ndet_ref, Ndet_non_ref,i_generator,n_selected,det_buffer,Nint,iproc,key_mask) @@ -75,11 +75,10 @@ subroutine mrcc_dress(delta_ij_, delta_ii_, Ndet_ref, Ndet_non_ref,i_generator,n integer :: i,j,k,l integer :: degree_alpha(psi_det_size) integer :: idx_alpha(0:psi_det_size) - logical :: good + logical :: good, fullMatch integer(bit_kind) :: tq(Nint,2,n_selected) integer :: N_tq, c_ref ,degree - integer :: connected_to_ref double precision :: hIk, hla, hIl, dIk(N_states), dka(N_states), dIa(N_states) double precision, allocatable :: dIa_hla(:,:) @@ -91,57 +90,20 @@ subroutine mrcc_dress(delta_ij_, delta_ii_, Ndet_ref, Ndet_non_ref,i_generator,n integer :: iint, ipos integer :: i_state, k_sd, l_sd, i_I, i_alpha - integer(bit_kind),allocatable :: miniList(:,:,:), supalist(:,:,:) + integer(bit_kind),allocatable :: miniList(:,:,:) integer(bit_kind),intent(in) :: key_mask(Nint, 2) integer,allocatable :: idx_miniList(:) - integer :: N_miniList, N_supalist, ni, leng + integer :: N_miniList, ni, leng leng = max(N_det_generators, N_det_non_ref) - allocate(miniList(Nint, 2, leng), idx_miniList(leng), supalist(Nint,2,leng)) + allocate(miniList(Nint, 2, leng), idx_miniList(leng)) - l = 0 - N_miniList = 0 - N_supalist = 0 - - do ni = 1,Nint - l += popcnt(key_mask(ni,1)) + popcnt(key_mask(ni,2)) - end do + !create_minilist_find_previous(key_mask, fullList, miniList, N_fullList, N_miniList, fullMatch, Nint) + call create_minilist_find_previous(key_mask, psi_det_generators, miniList, i_generator-1, N_miniList, fullMatch, Nint) - if(l == 0) then - N_miniList = i_generator-1 - miniList(:,:,:N_miniList) = psi_det_generators(:,:,:N_minilist) - else - do i=i_generator-1,1,-1 - k = l - do ni=1,nint - k -= popcnt(iand(key_mask(ni,1), psi_det_generators(ni,1,i))) + popcnt(iand(key_mask(ni,2), psi_det_generators(ni,2,i))) - end do - -! if(k == 0) then -! deallocate(miniList, supalist, idx_miniList) -! return -! else if(k <= 2) then -! N_minilist += 1 -! miniList(:,:,N_minilist) = psi_det_generators(:,:,i) -! end if -! - if(k == 2) then - N_supalist += 1 - supalist(:,:,N_supalist) = psi_det_generators(:,:,i) - else if(k == 1) then - N_minilist += 1 - miniList(:,:,N_minilist) = psi_det_generators(:,:,i) - else if(k == 0) then - deallocate(miniList, supalist, idx_miniList) - return - end if - end do - end if - - if(N_supalist > 0) then - miniList(:,:,N_minilist+1:N_minilist+N_supalist) = supalist(:,:,:N_supalist) - N_minilist = N_minilist + N_supalist + if(fullMatch) then + return end if @@ -299,6 +261,7 @@ subroutine find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq integer :: nt,ni + logical, external :: is_connected_to integer(bit_kind),intent(in) :: miniList(Nint,2,N_det_generators) @@ -310,15 +273,18 @@ subroutine find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq i_loop : do i=1,N_selected - do j=1,N_miniList - nt = 0 - do ni=1,Nint - nt += popcnt(xor(miniList(ni,1,j), det_buffer(ni,1,i))) + popcnt(xor(miniList(ni,2,j), det_buffer(ni,2,i))) - end do - if(nt <= 4) then - cycle i_loop - end if - end do + if(is_connected_to(det_buffer(ni,1,i), miniList, Nint, N_miniList)) then + cycle + end if +! do j=1,N_miniList +! nt = 0 +! do ni=1,Nint +! nt += popcnt(xor(miniList(ni,1,j), det_buffer(ni,1,i))) + popcnt(xor(miniList(ni,2,j), det_buffer(ni,2,i))) +! end do +! if(nt <= 4) then +! cycle i_loop +! end if +! end do ! if(connected_to_ref(det_buffer(1,1,i),psi_det_generators,Nint, & ! i_generator,N_det_generators) /= 0) then ! cycle i_loop diff --git a/plugins/Perturbation/Moller_plesset.irp.f b/plugins/Perturbation/Moller_plesset.irp.f index 7435f70c..2e8ba8e1 100644 --- a/plugins/Perturbation/Moller_plesset.irp.f +++ b/plugins/Perturbation/Moller_plesset.irp.f @@ -31,7 +31,8 @@ subroutine pt2_moller_plesset(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_s (Fock_matrix_diag_mo(p1) + Fock_matrix_diag_mo(p2)) delta_e = 1.d0/delta_e - call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det,psi_selectors_size,n_st,i_H_psi_array) + !call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det,psi_selectors_size,n_st,i_H_psi_array) + call i_H_psi_nominilist(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det,psi_selectors_size,n_st,i_H_psi_array) h = diag_H_mat_elem(det_pert,Nint) do i =1,n_st H_pert_diag(i) = h diff --git a/plugins/Perturbation/delta_rho_perturbation.irp.f b/plugins/Perturbation/delta_rho_perturbation.irp.f index d83eb9a8..77972c88 100644 --- a/plugins/Perturbation/delta_rho_perturbation.irp.f +++ b/plugins/Perturbation/delta_rho_perturbation.irp.f @@ -1,4 +1,4 @@ -subroutine pt2_delta_rho_one_point(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_st) +subroutine pt2_delta_rho_one_point(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_st,minilist,idx_minilist,N_minilist) use bitmasks implicit none integer, intent(in) :: Nint,ndet,n_st @@ -7,6 +7,10 @@ subroutine pt2_delta_rho_one_point(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,nde double precision :: i_O1_psi_array(N_st) double precision :: i_H_psi_array(N_st) + integer, intent(in) :: N_minilist + integer, intent(in) :: idx_minilist(0:N_det_selectors) + integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors) + BEGIN_DOC ! compute the perturbatibe contribution to the Integrated Spin density at z = z_one point of one determinant ! @@ -46,7 +50,8 @@ subroutine pt2_delta_rho_one_point(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,nde ! endif call i_O1_psi_alpha_beta(mo_integrated_delta_rho_one_point,det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_O1_psi_array) - call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array) + !call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array) + call i_H_psi(det_pert,minilist,idx_minilist,N_minilist,psi_selectors_coef,Nint,N_minilist,psi_selectors_size,N_st,i_H_psi_array) h = diag_H_mat_elem(det_pert,Nint) oii = diag_O1_mat_elem_alpha_beta(mo_integrated_delta_rho_one_point,det_pert,N_int) diff --git a/plugins/Perturbation/dipole_moment.irp.f b/plugins/Perturbation/dipole_moment.irp.f index ca09c31c..4af9ea6b 100644 --- a/plugins/Perturbation/dipole_moment.irp.f +++ b/plugins/Perturbation/dipole_moment.irp.f @@ -1,4 +1,4 @@ -subroutine pt2_dipole_moment_z(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_st) +subroutine pt2_dipole_moment_z(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_st,minilist,idx_minilist,N_minilist) use bitmasks implicit none integer, intent(in) :: Nint,ndet,n_st @@ -7,6 +7,10 @@ subroutine pt2_dipole_moment_z(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_ double precision :: i_O1_psi_array(N_st) double precision :: i_H_psi_array(N_st) + integer, intent(in) :: N_minilist + integer, intent(in) :: idx_minilist(0:N_det_selectors) + integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors) + BEGIN_DOC ! compute the perturbatibe contribution to the dipole moment of one determinant ! @@ -46,7 +50,9 @@ subroutine pt2_dipole_moment_z(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_ ! endif call i_O1_psi(mo_dipole_z,det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_O1_psi_array) - call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array) + !call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array) + call i_H_psi(det_pert,minilist,idx_minilist,N_minilist,psi_selectors_coef,Nint,N_minilist,psi_selectors_size,N_st,i_H_psi_array) + h = diag_H_mat_elem(det_pert,Nint) oii = diag_O1_mat_elem(mo_dipole_z,det_pert,N_int) diff --git a/plugins/Perturbation/epstein_nesbet.irp.f b/plugins/Perturbation/epstein_nesbet.irp.f index 62cb0cd6..b4ad9bfe 100644 --- a/plugins/Perturbation/epstein_nesbet.irp.f +++ b/plugins/Perturbation/epstein_nesbet.irp.f @@ -1,4 +1,4 @@ -subroutine pt2_epstein_nesbet(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st) +subroutine pt2_epstein_nesbet(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st,minilist,idx_minilist,N_minilist) use bitmasks implicit none integer, intent(in) :: Nint,ndet,N_st @@ -6,6 +6,10 @@ subroutine pt2_epstein_nesbet(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_s double precision , intent(out) :: c_pert(N_st),e_2_pert(N_st),H_pert_diag(N_st) double precision :: i_H_psi_array(N_st) + integer, intent(in) :: N_minilist + integer, intent(in) :: idx_minilist(0:N_det_selectors) + integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors) + BEGIN_DOC ! compute the standard Epstein-Nesbet perturbative first order coefficient and second order energetic contribution ! @@ -23,7 +27,10 @@ subroutine pt2_epstein_nesbet(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_s ASSERT (Nint == N_int) ASSERT (Nint > 0) - call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array) + !call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array) + call i_H_psi(det_pert,minilist,idx_minilist,N_minilist,psi_selectors_coef,Nint,N_minilist,psi_selectors_size,N_st,i_H_psi_array) + + h = diag_H_mat_elem(det_pert,Nint) do i =1,N_st if(CI_electronic_energy(i)>h.and.CI_electronic_energy(i).ne.0.d0)then @@ -42,7 +49,7 @@ subroutine pt2_epstein_nesbet(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_s end -subroutine pt2_epstein_nesbet_2x2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st) +subroutine pt2_epstein_nesbet_2x2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st,minilist,idx_minilist,N_minilist) use bitmasks implicit none integer, intent(in) :: Nint,ndet,N_st @@ -50,6 +57,10 @@ subroutine pt2_epstein_nesbet_2x2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet double precision , intent(out) :: c_pert(N_st),e_2_pert(N_st),H_pert_diag(N_st) double precision :: i_H_psi_array(N_st) + integer, intent(in) :: N_minilist + integer, intent(in) :: idx_minilist(0:N_det_selectors) + integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors) + BEGIN_DOC ! compute the Epstein-Nesbet 2x2 diagonalization coefficient and energetic contribution ! @@ -67,7 +78,9 @@ subroutine pt2_epstein_nesbet_2x2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet ASSERT (Nint > 0) PROVIDE CI_electronic_energy - call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array) + !call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array) + call i_H_psi(det_pert,minilist,idx_minilist,N_minilist,psi_selectors_coef,Nint,N_minilist,psi_selectors_size,N_st,i_H_psi_array) + h = diag_H_mat_elem(det_pert,Nint) do i =1,N_st if (i_H_psi_array(i) /= 0.d0) then diff --git a/plugins/Perturbation/pert_sc2.irp.f b/plugins/Perturbation/pert_sc2.irp.f index bdd8f97c..15399e4e 100644 --- a/plugins/Perturbation/pert_sc2.irp.f +++ b/plugins/Perturbation/pert_sc2.irp.f @@ -1,5 +1,5 @@ -subroutine pt2_epstein_nesbet_SC2_projected(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st) +subroutine pt2_epstein_nesbet_SC2_projected(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st,minilist,idx_minilist,N_minilist) use bitmasks implicit none integer, intent(in) :: Nint,ndet,N_st @@ -8,6 +8,10 @@ subroutine pt2_epstein_nesbet_SC2_projected(det_pert,c_pert,e_2_pert,H_pert_diag double precision :: i_H_psi_array(N_st) integer :: idx_repeat(0:ndet) + integer, intent(in) :: N_minilist + integer, intent(in) :: idx_minilist(0:N_det_selectors) + integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors) + BEGIN_DOC ! compute the Epstein-Nesbet perturbative first order coefficient and second order energetic contribution ! @@ -84,7 +88,7 @@ subroutine pt2_epstein_nesbet_SC2_projected(det_pert,c_pert,e_2_pert,H_pert_diag end -subroutine pt2_epstein_nesbet_SC2_no_projected(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st) +subroutine pt2_epstein_nesbet_SC2_no_projected(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st,minilist,idx_minilist,N_minilist) use bitmasks implicit none integer, intent(in) :: Nint,ndet,N_st @@ -93,6 +97,10 @@ subroutine pt2_epstein_nesbet_SC2_no_projected(det_pert,c_pert,e_2_pert,H_pert_d double precision :: i_H_psi_array(N_st) integer :: idx_repeat(0:ndet) + integer, intent(in) :: N_minilist + integer, intent(in) :: idx_minilist(0:N_det_selectors) + integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors) + BEGIN_DOC ! compute the Epstein-Nesbet perturbative first order coefficient and second order energetic contribution ! @@ -183,7 +191,7 @@ double precision function repeat_all_e_corr(key_in) end -subroutine pt2_epstein_nesbet_sc2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st) +subroutine pt2_epstein_nesbet_sc2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st,minilist,idx_minilist,N_minilist) use bitmasks implicit none integer, intent(in) :: Nint,ndet,N_st @@ -191,6 +199,10 @@ subroutine pt2_epstein_nesbet_sc2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet double precision , intent(out) :: c_pert(N_st),e_2_pert(N_st),H_pert_diag(N_st) double precision :: i_H_psi_array(N_st) + integer, intent(in) :: N_minilist + integer, intent(in) :: idx_minilist(0:N_det_selectors) + integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors) + BEGIN_DOC ! compute the standard Epstein-Nesbet perturbative first order coefficient and second order energetic contribution ! @@ -208,7 +220,10 @@ subroutine pt2_epstein_nesbet_sc2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet ASSERT (Nint == N_int) ASSERT (Nint > 0) - call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array) + !call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array) + call i_H_psi(det_pert,minilist,idx_minilist,N_minilist,psi_selectors_coef,Nint,N_minilist,psi_selectors_size,N_st,i_H_psi_array) + + h = diag_H_mat_elem(det_pert,Nint) do i =1,N_st if(CI_SC2_electronic_energy(i)>h.and.CI_SC2_electronic_energy(i).ne.0.d0)then diff --git a/plugins/Perturbation/pert_single.irp.f b/plugins/Perturbation/pert_single.irp.f index d04ca7ca..e2fbc9bf 100644 --- a/plugins/Perturbation/pert_single.irp.f +++ b/plugins/Perturbation/pert_single.irp.f @@ -1,4 +1,4 @@ -subroutine pt2_h_core(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st) +subroutine pt2_h_core(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st,minilist,idx_minilist,N_minilist) use bitmasks implicit none integer, intent(in) :: Nint,ndet,N_st @@ -6,6 +6,10 @@ subroutine pt2_h_core(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st) double precision , intent(out) :: c_pert(N_st),e_2_pert(N_st),H_pert_diag(N_st) double precision :: i_H_psi_array(N_st) + integer, intent(in) :: N_minilist + integer, intent(in) :: idx_minilist(0:N_det_selectors) + integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors) + BEGIN_DOC ! compute the standard Epstein-Nesbet perturbative first order coefficient and second order energetic contribution ! diff --git a/plugins/Perturbation/perturbation.template.f b/plugins/Perturbation/perturbation.template.f index a5ab12e7..05176fe6 100644 --- a/plugins/Perturbation/perturbation.template.f +++ b/plugins/Perturbation/perturbation.template.f @@ -2,7 +2,7 @@ BEGIN_SHELL [ /usr/bin/env python ] import perturbation END_SHELL -subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,coef_pert_buffer,sum_e_2_pert,sum_norm_pert,sum_H_pert_diag,N_st,Nint) +subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,coef_pert_buffer,sum_e_2_pert,sum_norm_pert,sum_H_pert_diag,N_st,Nint,key_mask) implicit none BEGIN_DOC ! Applly pertubration ``$PERT`` to the buffer of determinants generated in the H_apply @@ -11,25 +11,59 @@ subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,c integer, intent(in) :: Nint, N_st, buffer_size, i_generator integer(bit_kind), intent(in) :: buffer(Nint,2,buffer_size) + integer(bit_kind),intent(in) :: key_mask(Nint,2) double precision, intent(inout) :: sum_norm_pert(N_st),sum_e_2_pert(N_st) 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 :: i,k, c_ref, ni, ex integer, external :: connected_to_ref logical, external :: is_in_wavefunction + integer(bit_kind) :: minilist(Nint,2,N_det_selectors) + integer :: idx_minilist(N_det_selectors), N_minilist + + integer(bit_kind) :: minilist_gen(Nint,2,N_det_generators) + integer :: N_minilist_gen + logical :: fullMatch + logical, external :: is_connected_to + + + ASSERT (Nint > 0) ASSERT (Nint == N_int) ASSERT (buffer_size >= 0) ASSERT (minval(sum_norm_pert) >= 0.d0) ASSERT (N_st > 0) - do i = 1,buffer_size + + call create_minilist(key_mask, psi_selectors, miniList, idx_miniList, N_det_selectors, N_minilist, Nint) + call create_minilist_find_previous(key_mask, psi_det_generators, miniList_gen, i_generator-1, N_minilist_gen, fullMatch, Nint) + + if(fullMatch) then + return + end if + + + buffer_loop : do i = 1,buffer_size + +! do k=1,N_minilist_gen +! ex = 0 +! do ni=1,Nint +! ex += popcnt(xor(minilist_gen(ni,1,k), buffer(ni,1,i))) + popcnt(xor(minilist_gen(ni,2,k), buffer(ni,2,i))) +! end do +! if(ex <= 4) then +! cycle buffer_loop +! end if +! end do + +! c_ref = connected_to_ref(buffer(1,1,i),miniList_gen,Nint,N_minilist_gen+1,N_minilist_gen) +! +! if (c_ref /= 0) then +! cycle +! endif - c_ref = connected_to_ref(buffer(1,1,i),psi_det_generators,Nint,i_generator,N_det_generators) - - if (c_ref /= 0) then + if(is_connected_to(buffer(1,1,i), miniList_gen, Nint, N_minilist_gen)) then cycle - endif + end if if (is_in_wavefunction(buffer(1,1,i),Nint)) then cycle @@ -37,8 +71,10 @@ subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,c integer :: degree call get_excitation_degree(HF_bitmask,buffer(1,1,i),degree,N_int) +! call pt2_$PERT(buffer(1,1,i), & +! c_pert,e_2_pert,H_pert_diag,Nint,N_det_selectors,n_st,minilist,idx_minilist) call pt2_$PERT(buffer(1,1,i), & - c_pert,e_2_pert,H_pert_diag,Nint,N_det_selectors,n_st) + c_pert,e_2_pert,H_pert_diag,Nint,N_minilist,n_st,minilist,idx_minilist,N_minilist) !!!!!!!!!!!!!!!!! MAUVAISE SIGNATURE PR LES AUTRES PT2_* !!!!! do k = 1,N_st e_2_pert_buffer(k,i) = e_2_pert(k) @@ -48,7 +84,7 @@ subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,c sum_H_pert_diag(k) += H_pert_diag(k) enddo - enddo + enddo buffer_loop end diff --git a/scripts/generate_h_apply.py b/scripts/generate_h_apply.py index b22797f9..4e95e0f1 100755 --- a/scripts/generate_h_apply.py +++ b/scripts/generate_h_apply.py @@ -205,7 +205,7 @@ class H_apply(object): """ 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, & - sum_norm_pert,sum_H_pert_diag,N_st,N_int) + sum_norm_pert,sum_H_pert_diag,N_st,N_int,key_mask) """%(pert,) self.data["finalization"] = """ """ diff --git a/src/Determinants/connected_to_ref.irp.f b/src/Determinants/connected_to_ref.irp.f index 8f594738..dc7698b5 100644 --- a/src/Determinants/connected_to_ref.irp.f +++ b/src/Determinants/connected_to_ref.irp.f @@ -154,6 +154,41 @@ integer function get_index_in_psi_det_sorted_bit(key,Nint) ! END DEBUG is_in_wf end + +logical function is_connected_to(key,keys,Nint,Ndet) + use bitmasks + implicit none + integer, intent(in) :: Nint, Ndet + integer(bit_kind), intent(in) :: keys(Nint,2,Ndet) + integer(bit_kind), intent(in) :: key(Nint,2) + + integer :: i, l + integer :: degree_x2 + + + ASSERT (Nint > 0) + ASSERT (Nint == N_int) + + is_connected_to = .false. + + do i=1,Ndet + degree_x2 = popcnt(xor( key(1,1), keys(1,1,i))) + & + popcnt(xor( key(1,2), keys(1,2,i))) + !DEC$ LOOP COUNT MIN(3) + do l=2,Nint + degree_x2 = degree_x2 + popcnt(xor( key(l,1), keys(l,1,i))) +& + popcnt(xor( key(l,2), keys(l,2,i))) + enddo + if (degree_x2 > 4) then + cycle + else + is_connected_to = .true. + return + endif + enddo +end + + integer function connected_to_ref(key,keys,Nint,N_past_in,Ndet) use bitmasks implicit none diff --git a/src/Determinants/slater_rules.irp.f b/src/Determinants/slater_rules.irp.f index 945450e2..954d2f6a 100644 --- a/src/Determinants/slater_rules.irp.f +++ b/src/Determinants/slater_rules.irp.f @@ -763,8 +763,107 @@ subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble) end +subroutine create_minilist(key_mask, fullList, miniList, idx_miniList, N_fullList, N_miniList, Nint) + use bitmasks + implicit none + + integer(bit_kind), intent(in) :: fullList(Nint, 2, N_fullList) + integer, intent(in) :: N_fullList + integer(bit_kind),intent(out) :: miniList(Nint, 2, N_fullList) + integer,intent(out) :: idx_miniList(N_fullList), N_miniList + integer, intent(in) :: Nint + integer(bit_kind) :: key_mask(Nint, 2) + integer :: ni, i, n_a, n_b, e_a, e_b + + + n_a = 0 + n_b = 0 + do ni=1,nint + n_a = n_a + popcnt(key_mask(ni,1)) + n_b = n_b + popcnt(key_mask(ni,2)) + end do + + if(n_a == 0) then + N_miniList = N_fullList + miniList(:,:,:) = fullList(:,:,:) + do i=1,N_fullList + idx_miniList(i) = i + end do + return + end if + + N_miniList = 0 + + do i=1,N_fullList + e_a = n_a + e_b = n_b + do ni=1,nint + e_a -= popcnt(iand(fullList(ni, 1, i), key_mask(ni, 1))) + e_b -= popcnt(iand(fullList(ni, 2, i), key_mask(ni, 2))) + end do + + if(e_a + e_b <= 2) then + N_miniList = N_miniList + 1 + miniList(:,:,N_miniList) = fullList(:,:,i) + idx_miniList(N_miniList) = i + end if + end do +end subroutine -subroutine i_H_psi(key,keys,coef,Nint,Ndet,Ndet_max,Nstate,i_H_psi_array) +subroutine create_minilist_find_previous(key_mask, fullList, miniList, N_fullList, N_miniList, fullMatch, Nint) + use bitmasks + implicit none + + integer(bit_kind), intent(in) :: fullList(Nint, 2, N_fullList) + integer, intent(in) :: N_fullList + integer(bit_kind),intent(out) :: miniList(Nint, 2, N_fullList) + integer(bit_kind) :: subList(Nint, 2, N_fullList) + logical,intent(out) :: fullMatch + integer,intent(out) :: N_miniList + integer, intent(in) :: Nint + integer(bit_kind) :: key_mask(Nint, 2) + integer :: ni, i, k, l, N_subList + + + fullMatch = .false. + l = 0 + N_miniList = 0 + N_subList = 0 + + do ni = 1,Nint + l += popcnt(key_mask(ni,1)) + popcnt(key_mask(ni,2)) + end do + + if(l == 0) then + N_miniList = N_fullList + miniList(:,:,:N_miniList) = fullList(:,:,:N_minilist) + else + do i=N_fullList,1,-1 + k = l + do ni=1,nint + k -= popcnt(iand(key_mask(ni,1), fullList(ni,1,i))) + popcnt(iand(key_mask(ni,2), fullList(ni,2,i))) + end do + if(k == 2) then + N_subList += 1 + subList(:,:,N_subList) = fullList(:,:,i) + else if(k == 1) then + N_minilist += 1 + miniList(:,:,N_minilist) = fullList(:,:,i) + else if(k == 0) then + fullMatch = .true. + return + end if + end do + end if + + if(N_subList > 0) then + miniList(:,:,N_minilist+1:N_minilist+N_subList) = sublist(:,:,:N_subList) + N_minilist = N_minilist + N_subList + end if +end subroutine + + +subroutine i_H_psi_nominilist(key,keys,coef,Nint,Ndet,Ndet_max,Nstate,i_H_psi_array) use bitmasks implicit none integer, intent(in) :: Nint, Ndet,Ndet_max,Nstate @@ -800,6 +899,50 @@ subroutine i_H_psi(key,keys,coef,Nint,Ndet,Ndet_max,Nstate,i_H_psi_array) enddo end + +subroutine i_H_psi(key,keys,idx_key,N_minilist,coef,Nint,Ndet,Ndet_max,Nstate,i_H_psi_array) + use bitmasks + implicit none + integer, intent(in) :: Nint, Ndet,Ndet_max,Nstate,idx_key(Ndet), N_minilist + integer(bit_kind), intent(in) :: keys(Nint,2,Ndet) + integer(bit_kind), intent(in) :: key(Nint,2) + double precision, intent(in) :: coef(Ndet_max,Nstate) + double precision, intent(out) :: i_H_psi_array(Nstate) + + integer :: i, ii,j, i_in_key, i_in_coef + double precision :: phase + integer :: exc(0:2,2,2) + double precision :: hij + integer :: idx(0:Ndet) + BEGIN_DOC + ! for the various Nstates + END_DOC + + ASSERT (Nint > 0) + ASSERT (N_int == Nint) + ASSERT (Nstate > 0) + ASSERT (Ndet > 0) + ASSERT (Ndet_max >= Ndet) + i_H_psi_array = 0.d0 + + !call filter_connected_i_H_psi0(keys,key,Nint,Ndet,idx) + call filter_connected_i_H_psi0(keys,key,Nint,N_minilist,idx) + do ii=1,idx(0) + !i = idx_key(idx(ii)) + i_in_key = idx(ii) + i_in_coef = idx_key(idx(ii)) + !DEC$ FORCEINLINE +! ! call i_H_j(keys(1,1,i),key,Nint,hij) +! ! do j = 1, Nstate +! ! i_H_psi_array(j) = i_H_psi_array(j) + coef(i,j)*hij +! ! enddo + call i_H_j(keys(1,1,i_in_key),key,Nint,hij) + do j = 1, Nstate + i_H_psi_array(j) = i_H_psi_array(j) + coef(i_in_coef,j)*hij + enddo + enddo +end + subroutine i_H_psi_sec_ord(key,keys,coef,Nint,Ndet,Ndet_max,Nstate,i_H_psi_array,idx_interaction,interactions) use bitmasks implicit none