From 007e234997f568a38be3dd8983c62c908663dd86 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Fri, 20 Nov 2015 19:51:56 +0100 Subject: [PATCH] bitstring_to_list_ab --- plugins/MP2/mp2_wf.irp.f | 2 +- scripts/generate_h_apply.py | 2 +- src/Determinants/density_matrix.irp.f | 2 - src/Determinants/slater_rules.irp.f | 88 +++++++++++++++++++-------- 4 files changed, 63 insertions(+), 31 deletions(-) diff --git a/plugins/MP2/mp2_wf.irp.f b/plugins/MP2/mp2_wf.irp.f index f751209b..ad068b8a 100644 --- a/plugins/MP2/mp2_wf.irp.f +++ b/plugins/MP2/mp2_wf.irp.f @@ -18,13 +18,13 @@ program mp2_wf call H_apply_mp2_selection(pt2, norm_pert, H_pert_diag, N_st) psi_det = psi_det_sorted psi_coef = psi_coef_sorted + touch N_det psi_det psi_coef print*,'N_det = ',N_det print*,'-----' print *, 'PT2 = ', pt2(1) print *, 'E = ', HF_energy print *, 'E_before +PT2 = ', HF_energy+pt2(1) N_det = min(N_det,N_det_max) - touch N_det psi_det psi_coef call save_wavefunction call ezfio_set_mp2_energy(HF_energy+pt2(1)) deallocate(pt2,norm_pert,H_pert_diag) diff --git a/scripts/generate_h_apply.py b/scripts/generate_h_apply.py index d0f67550..26f19d0d 100755 --- a/scripts/generate_h_apply.py +++ b/scripts/generate_h_apply.py @@ -60,7 +60,7 @@ class H_apply(object): s["omp_master"] = "!$OMP MASTER" s["omp_end_master"] = "!$OMP END MASTER" s["omp_barrier"] = "!$OMP BARRIER" - s["omp_do"] = "!$OMP DO SCHEDULE (static)" + s["omp_do"] = "!$OMP DO SCHEDULE (static,1)" s["omp_enddo"] = "!$OMP ENDDO NOWAIT" d = { True : '.True.', False : '.False.'} diff --git a/src/Determinants/density_matrix.irp.f b/src/Determinants/density_matrix.irp.f index f72b337c..5f087642 100644 --- a/src/Determinants/density_matrix.irp.f +++ b/src/Determinants/density_matrix.irp.f @@ -71,7 +71,6 @@ one_body_dm_mo_beta = one_body_dm_mo_beta + tmp_b !$OMP END CRITICAL deallocate(tmp_a,tmp_b) - !$OMP BARRIER !$OMP END PARALLEL endif @@ -157,7 +156,6 @@ END_PROVIDER one_body_single_double_dm_mo_beta = one_body_single_double_dm_mo_beta + tmp_b !$OMP END CRITICAL deallocate(tmp_a,tmp_b) - !$OMP BARRIER !$OMP END PARALLEL END_PROVIDER diff --git a/src/Determinants/slater_rules.irp.f b/src/Determinants/slater_rules.irp.f index b60bf0f7..5337670f 100644 --- a/src/Determinants/slater_rules.irp.f +++ b/src/Determinants/slater_rules.irp.f @@ -349,6 +349,42 @@ subroutine get_mono_excitation(det1,det2,exc,phase,Nint) enddo end +subroutine bitstring_to_list_ab( string, list, n_elements, Nint) + use bitmasks + implicit none + BEGIN_DOC + ! Gives the inidices(+1) of the bits set to 1 in the bit string + ! For alpha/beta determinants + END_DOC + integer, intent(in) :: Nint + integer(bit_kind), intent(in) :: string(Nint,2) + integer, intent(out) :: list(Nint*bit_kind_size,2) + integer, intent(out) :: n_elements(2) + + integer :: i, ishift + integer(bit_kind) :: l + + n_elements(1) = 0 + n_elements(2) = 0 + ishift = 2 + do i=1,Nint + l = string(i,1) + do while (l /= 0_bit_kind) + n_elements(1) = n_elements(1)+1 + list(n_elements(1),1) = ishift+popcnt(l-1_bit_kind) - popcnt(l) + l = iand(l,l-1_bit_kind) + enddo + l = string(i,2) + do while (l /= 0_bit_kind) + n_elements(2) = n_elements(2)+1 + list(n_elements(2),2) = ishift+popcnt(l-1_bit_kind) - popcnt(l) + l = iand(l,l-1_bit_kind) + enddo + ishift = ishift + bit_kind_size + enddo + +end + @@ -370,7 +406,7 @@ subroutine i_H_j(key_i,key_j,Nint,hij) integer :: i,j,k integer :: occ(Nint*bit_kind_size,2) double precision :: diag_H_mat_elem, phase,phase_2 - integer :: n_occ_alpha, n_occ_beta + integer :: n_occ_ab(2) logical :: has_mipi(Nint*bit_kind_size) double precision :: mipi(Nint*bit_kind_size), miip(Nint*bit_kind_size) PROVIDE mo_bielec_integrals_in_map mo_integrals_map @@ -422,8 +458,8 @@ subroutine i_H_j(key_i,key_j,Nint,hij) endif case (1) call get_mono_excitation(key_i,key_j,exc,phase,Nint) - call bitstring_to_list(key_i(1,1), occ(1,1), n_occ_alpha, Nint) - call bitstring_to_list(key_i(1,2), occ(1,2), n_occ_beta, Nint) + !DIR$ FORCEINLINE + call bitstring_to_list_ab(key_i, occ, n_occ_ab, Nint) has_mipi = .False. if (exc(0,1,1) == 1) then ! Mono alpha @@ -506,7 +542,7 @@ subroutine i_H_j_phase_out(key_i,key_j,Nint,hij,phase,exc,degree) integer :: i,j,k integer :: occ(Nint*bit_kind_size,2) double precision :: diag_H_mat_elem - integer :: n_occ_alpha, n_occ_beta + integer :: n_occ_ab(2) logical :: has_mipi(Nint*bit_kind_size) double precision :: mipi(Nint*bit_kind_size), miip(Nint*bit_kind_size) PROVIDE mo_bielec_integrals_in_map mo_integrals_map @@ -558,8 +594,8 @@ subroutine i_H_j_phase_out(key_i,key_j,Nint,hij,phase,exc,degree) endif case (1) call get_mono_excitation(key_i,key_j,exc,phase,Nint) - call bitstring_to_list(key_i(1,1), occ(1,1), n_occ_alpha, Nint) - call bitstring_to_list(key_i(1,2), occ(1,2), n_occ_beta, Nint) + !DIR$ FORCEINLINE + call bitstring_to_list_ab(key_i, occ, n_occ_ab, Nint) has_mipi = .False. if (exc(0,1,1) == 1) then ! Mono alpha @@ -642,7 +678,7 @@ subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble) integer :: i,j,k integer :: occ(Nint*bit_kind_size,2) double precision :: diag_H_mat_elem, phase,phase_2 - integer :: n_occ_alpha, n_occ_beta + integer :: n_occ_ab(2) logical :: has_mipi(Nint*bit_kind_size) double precision :: mipi(Nint*bit_kind_size), miip(Nint*bit_kind_size) PROVIDE mo_bielec_integrals_in_map mo_integrals_map @@ -696,8 +732,8 @@ subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble) endif case (1) call get_mono_excitation(key_i,key_j,exc,phase,Nint) - call bitstring_to_list(key_i(1,1), occ(1,1), n_occ_alpha, Nint) - call bitstring_to_list(key_i(1,2), occ(1,2), n_occ_beta, Nint) + !DIR$ FORCEINLINE + call bitstring_to_list_ab(key_i, occ, n_occ_ab, Nint) has_mipi = .False. if (exc(0,1,1) == 1) then ! Mono alpha @@ -1229,15 +1265,15 @@ double precision function diag_H_mat_elem(det_in,Nint) endif !call debug_det(det_in,Nint) - integer :: tmp - call bitstring_to_list(particle(1,1), occ_particle(1,1), tmp, Nint) - ASSERT (tmp == nexc(1)) - call bitstring_to_list(particle(1,2), occ_particle(1,2), tmp, Nint) - ASSERT (tmp == nexc(2)) - call bitstring_to_list(hole(1,1), occ_hole(1,1), tmp, Nint) - ASSERT (tmp == nexc(1)) - call bitstring_to_list(hole(1,2), occ_hole(1,2), tmp, Nint) - ASSERT (tmp == nexc(2)) + integer :: tmp(2) + !DIR$ FORCEINLINE + call bitstring_to_list_ab(particle, occ_particle, tmp, Nint) + ASSERT (tmp(1) == nexc(1)) + ASSERT (tmp(2) == nexc(2)) + !DIR$ FORCEINLINE + call bitstring_to_list_ab(hole, occ_hole, tmp, Nint) + ASSERT (tmp(1) == nexc(1)) + ASSERT (tmp(2) == nexc(2)) det_tmp = ref_bitmask do ispin=1,2 @@ -1317,13 +1353,11 @@ subroutine ac_operator(iorb,ispin,key,hjj,Nint,na,nb) ASSERT (ispin < 3) ASSERT (Nint > 0) - integer :: tmp + integer :: tmp(2) !DIR$ FORCEINLINE - call bitstring_to_list(key(1,1), occ(1,1), tmp, Nint) - ASSERT (tmp == elec_alpha_num) - !DIR$ FORCEINLINE - call bitstring_to_list(key(1,2), occ(1,2), tmp, Nint) - ASSERT (tmp == elec_beta_num) + call bitstring_to_list_ab(key, occ, tmp, Nint) + ASSERT (tmp(1) == elec_alpha_num) + ASSERT (tmp(2) == elec_beta_num) k = ishft(iorb-1,-bit_kind_shift)+1 ASSERT (k > 0) @@ -1354,10 +1388,10 @@ subroutine get_occ_from_key(key,occ,Nint) integer(bit_kind), intent(in) :: key(Nint,2) integer , intent(in) :: Nint integer , intent(out) :: occ(Nint*bit_kind_size,2) - integer :: tmp + integer :: tmp(2) - call bitstring_to_list(key(1,1), occ(1,1), tmp, Nint) - call bitstring_to_list(key(1,2), occ(1,2), tmp, Nint) + !DIR$ FORCEINLINE + call bitstring_to_list_ab(key, occ, tmp, Nint) end