From d13853691a2038f61c3d13b9931b82af96faa50d Mon Sep 17 00:00:00 2001 From: Emmanuel Giner Date: Wed, 2 Nov 2016 17:39:39 +0100 Subject: [PATCH] conflicts minimized for merge --- config/ifort.cfg | 9 +- plugins/Full_CI/EZFIO.cfg | 12 ++ plugins/Full_CI/full_ci.irp.f | 5 +- plugins/Full_CI/full_ci_no_skip.irp.f | 4 + plugins/MRPT_Utils/excitations_cas.irp.f | 40 +++---- plugins/MRPT_Utils/fock_like_operators.irp.f | 12 +- plugins/MRPT_Utils/new_way.irp.f | 26 ++--- .../new_way_second_order_coef.irp.f | 12 +- plugins/MRPT_Utils/second_order_new.irp.f | 12 +- plugins/MRPT_Utils/second_order_new_2p.irp.f | 6 +- plugins/Perturbation/EZFIO.cfg | 1 + plugins/Properties/EZFIO.cfg | 1 + plugins/Properties/hyperfine_constants.irp.f | 2 +- plugins/Selectors_full/selectors.irp.f | 3 +- src/Bitmask/bitmask_cas_routines.irp.f | 24 ++-- src/Determinants/EZFIO.cfg | 4 +- src/Determinants/NEEDED_CHILDREN_MODULES | 2 +- src/Determinants/diagonalize_CI.irp.f | 12 -- src/Determinants/slater_rules.irp.f | 104 +++++++++--------- src/Determinants/test_3d.irp.f | 15 --- src/Determinants/test_two_body.irp.f | 18 --- src/Determinants/utils.irp.f | 35 +----- src/Integrals_Bielec/integrals_3_index.irp.f | 6 +- src/Integrals_Bielec/map_integrals.irp.f | 26 +---- src/Integrals_Bielec/mo_bi_integrals.irp.f | 20 +--- 25 files changed, 155 insertions(+), 256 deletions(-) delete mode 100644 src/Determinants/test_3d.irp.f delete mode 100644 src/Determinants/test_two_body.irp.f diff --git a/config/ifort.cfg b/config/ifort.cfg index da414912..c1d7e968 100644 --- a/config/ifort.cfg +++ b/config/ifort.cfg @@ -31,14 +31,14 @@ OPENMP : 1 ; Append OpenMP flags # -ftz : Flushes denormal results to zero # [OPT] -FCFLAGS : -O2 -xHost -ip -ftz +FCFLAGS : -xSSE4.2 -O2 -ip -opt-prefetch -ftz -g # Profiling flags ################# # [PROFILE] FC : -p -g -FCFLAGS : -xSSE4.2 -O2 -ip -ftz +FCFLAGS : -xSSE4.2 -O2 -ip -opt-prefetch -ftz # Debugging flags ################# @@ -50,9 +50,8 @@ FCFLAGS : -xSSE4.2 -O2 -ip -ftz # -xSSE2 : Valgrind needs a very simple x86 executable # [DEBUG] -FC : -g -traceback -fpe0 -FCFLAGS : -xSSE2 -C -IRPF90_FLAGS : --openmp +FC : -g -traceback +FCFLAGS : -xSSE2 -C -fpe0 # OpenMP flags ################# diff --git a/plugins/Full_CI/EZFIO.cfg b/plugins/Full_CI/EZFIO.cfg index 9a552cd0..afb25d2e 100644 --- a/plugins/Full_CI/EZFIO.cfg +++ b/plugins/Full_CI/EZFIO.cfg @@ -8,3 +8,15 @@ type: double precision doc: Calculated FCI energy + PT2 interface: ezfio +[threshold_generators_pt2] +type: Threshold +doc: Thresholds on generators (fraction of the norm) for final PT2 calculation +interface: ezfio,provider,ocaml +default: 0.999 + +[threshold_selectors_pt2] +type: Threshold +doc: Thresholds on selectors (fraction of the norm) for final PT2 calculation +interface: ezfio,provider,ocaml +default: 1. + diff --git a/plugins/Full_CI/full_ci.irp.f b/plugins/Full_CI/full_ci.irp.f index e6d0f7f2..ff599870 100644 --- a/plugins/Full_CI/full_ci.irp.f +++ b/plugins/Full_CI/full_ci.irp.f @@ -90,8 +90,9 @@ program full_ci call diagonalize_CI if(do_pt2_end)then print*,'Last iteration only to compute the PT2' - threshold_selectors = 1.d0 - threshold_generators = 0.999d0 + threshold_generators = threshold_generators_pt2 + threshold_selectors = threshold_selectors_pt2 + SOFT_TOUCH threshold_generators threshold_selectors call H_apply_FCI_PT2(pt2, norm_pert, H_pert_diag, N_st) print *, 'Final step' diff --git a/plugins/Full_CI/full_ci_no_skip.irp.f b/plugins/Full_CI/full_ci_no_skip.irp.f index 3ed304a1..078334f7 100644 --- a/plugins/Full_CI/full_ci_no_skip.irp.f +++ b/plugins/Full_CI/full_ci_no_skip.irp.f @@ -73,6 +73,10 @@ program full_ci call diagonalize_CI if(do_pt2_end)then print*,'Last iteration only to compute the PT2' + threshold_generators = threshold_generators_pt2 + threshold_selectors = threshold_selectors_pt2 + SOFT_TOUCH threshold_generators threshold_selectors + ! print*,'The thres' call H_apply_FCI_PT2(pt2, norm_pert, H_pert_diag, N_st) diff --git a/plugins/MRPT_Utils/excitations_cas.irp.f b/plugins/MRPT_Utils/excitations_cas.irp.f index 805070f7..fb5cc953 100644 --- a/plugins/MRPT_Utils/excitations_cas.irp.f +++ b/plugins/MRPT_Utils/excitations_cas.irp.f @@ -270,7 +270,7 @@ subroutine i_H_j_dyall(key_i,key_j,Nint,hij) integer :: exc(0:2,2,2) integer :: degree - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral integer :: m,n,p,q integer :: i,j,k integer :: occ(Nint*bit_kind_size,2) @@ -291,31 +291,31 @@ subroutine i_H_j_dyall(key_i,key_j,Nint,hij) call get_double_excitation(key_i,key_j,exc,phase,Nint) if (exc(0,1,1) == 1) then ! Mono alpha, mono beta - hij = phase*get_mo_bielec_integral_schwartz( & + hij = phase*get_mo_bielec_integral( & exc(1,1,1), & exc(1,1,2), & exc(1,2,1), & exc(1,2,2) ,mo_integrals_map) else if (exc(0,1,1) == 2) then ! Double alpha - hij = phase*(get_mo_bielec_integral_schwartz( & + hij = phase*(get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(1,2,1), & exc(2,2,1) ,mo_integrals_map) - & - get_mo_bielec_integral_schwartz( & + get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(2,2,1), & exc(1,2,1) ,mo_integrals_map) ) else if (exc(0,1,2) == 2) then ! Double beta - hij = phase*(get_mo_bielec_integral_schwartz( & + hij = phase*(get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(1,2,2), & exc(2,2,2) ,mo_integrals_map) - & - get_mo_bielec_integral_schwartz( & + get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(2,2,2), & @@ -333,15 +333,15 @@ subroutine i_H_j_dyall(key_i,key_j,Nint,hij) do k = 1, n_occ_ab(1) i = occ(k,1) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) - miip(i) = get_mo_bielec_integral_schwartz(m,i,i,p,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) + miip(i) = get_mo_bielec_integral(m,i,i,p,mo_integrals_map) has_mipi(i) = .True. endif enddo do k = 1, n_occ_ab(2) i = occ(k,2) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) has_mipi(i) = .True. endif enddo @@ -360,15 +360,15 @@ subroutine i_H_j_dyall(key_i,key_j,Nint,hij) do k = 1, n_occ_ab(2) i = occ(k,2) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) - miip(i) = get_mo_bielec_integral_schwartz(m,i,i,p,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) + miip(i) = get_mo_bielec_integral(m,i,i,p,mo_integrals_map) has_mipi(i) = .True. endif enddo do k = 1, n_occ_ab(1) i = occ(k,1) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) has_mipi(i) = .True. endif enddo @@ -494,7 +494,7 @@ subroutine i_H_j_dyall_no_exchange(key_i,key_j,Nint,hij) integer :: exc(0:2,2,2) integer :: degree - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral integer :: m,n,p,q integer :: i,j,k integer :: occ(Nint*bit_kind_size,2) @@ -518,7 +518,7 @@ subroutine i_H_j_dyall_no_exchange(key_i,key_j,Nint,hij) if(exc(1,1,1) == exc(1,2,2) .and. exc(1,2,1) == exc(1,1,2))then hij = 0.d0 else - hij = phase*get_mo_bielec_integral_schwartz( & + hij = phase*get_mo_bielec_integral( & exc(1,1,1), & exc(1,1,2), & exc(1,2,1), & @@ -526,14 +526,14 @@ subroutine i_H_j_dyall_no_exchange(key_i,key_j,Nint,hij) endif else if (exc(0,1,1) == 2) then ! Double alpha - hij = phase*get_mo_bielec_integral_schwartz( & + hij = phase*get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(1,2,1), & exc(2,2,1) ,mo_integrals_map) else if (exc(0,1,2) == 2) then ! Double beta - hij = phase*get_mo_bielec_integral_schwartz( & + hij = phase*get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(1,2,2), & @@ -551,14 +551,14 @@ subroutine i_H_j_dyall_no_exchange(key_i,key_j,Nint,hij) do k = 1, n_occ_ab(1) i = occ(k,1) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) has_mipi(i) = .True. endif enddo do k = 1, n_occ_ab(2) i = occ(k,2) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) has_mipi(i) = .True. endif enddo @@ -577,14 +577,14 @@ subroutine i_H_j_dyall_no_exchange(key_i,key_j,Nint,hij) do k = 1, n_occ_ab(2) i = occ(k,2) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) has_mipi(i) = .True. endif enddo do k = 1, n_occ_ab(1) i = occ(k,1) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) has_mipi(i) = .True. endif enddo diff --git a/plugins/MRPT_Utils/fock_like_operators.irp.f b/plugins/MRPT_Utils/fock_like_operators.irp.f index 2074daf6..d4ce0661 100644 --- a/plugins/MRPT_Utils/fock_like_operators.irp.f +++ b/plugins/MRPT_Utils/fock_like_operators.irp.f @@ -51,7 +51,7 @@ double precision :: accu_coulomb,accu_exchange(2) double precision :: na,nb,ntot double precision :: coulomb, exchange - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral integer :: j_act_orb,k_act_orb,i_inact_core_orb integer :: i_state @@ -75,8 +75,8 @@ na = one_body_dm_mo_alpha(j_act_orb,k_act_orb,i_state) nb = one_body_dm_mo_beta(j_act_orb,k_act_orb,i_state) ntot = na + nb - coulomb = get_mo_bielec_integral_schwartz(j_act_orb,i_inact_core_orb,k_act_orb,i_inact_core_orb,mo_integrals_map) - exchange = get_mo_bielec_integral_schwartz(j_act_orb,k_act_orb,i_inact_core_orb,i_inact_core_orb,mo_integrals_map) + coulomb = get_mo_bielec_integral(j_act_orb,i_inact_core_orb,k_act_orb,i_inact_core_orb,mo_integrals_map) + exchange = get_mo_bielec_integral(j_act_orb,k_act_orb,i_inact_core_orb,i_inact_core_orb,mo_integrals_map) accu_coulomb += 2.d0 * ntot * coulomb accu_exchange(1) += 2.d0 * na * exchange accu_exchange(2) += 2.d0 * nb * exchange @@ -97,7 +97,7 @@ double precision :: accu_coulomb,accu_exchange(2) double precision :: na,nb,ntot double precision :: coulomb, exchange - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral integer :: j_act_orb,i_virt_orb,k_act_orb integer :: i_state ! TODO : inverse loop of i_state @@ -122,8 +122,8 @@ na = one_body_dm_mo_alpha(j_act_orb,k_act_orb,i_state) nb = one_body_dm_mo_beta(j_act_orb,k_act_orb,i_state) ntot = na + nb - coulomb = get_mo_bielec_integral_schwartz(j_act_orb,i_virt_orb,k_act_orb,i_virt_orb,mo_integrals_map) - exchange = get_mo_bielec_integral_schwartz(j_act_orb,k_act_orb,i_virt_orb,i_virt_orb,mo_integrals_map) + coulomb = get_mo_bielec_integral(j_act_orb,i_virt_orb,k_act_orb,i_virt_orb,mo_integrals_map) + exchange = get_mo_bielec_integral(j_act_orb,k_act_orb,i_virt_orb,i_virt_orb,mo_integrals_map) accu_coulomb += 2.d0 * ntot * coulomb accu_exchange(1) += 2.d0 * na * exchange accu_exchange(2) += 2.d0 * nb * exchange diff --git a/plugins/MRPT_Utils/new_way.irp.f b/plugins/MRPT_Utils/new_way.irp.f index 09016ab0..fa5812e1 100644 --- a/plugins/MRPT_Utils/new_way.irp.f +++ b/plugins/MRPT_Utils/new_way.irp.f @@ -15,7 +15,7 @@ subroutine give_2h1p_contrib(matrix_2h1p) integer(bit_kind) :: det_tmp(N_int,2) integer :: exc(0:2,2,2) integer :: accu_elec - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral double precision :: active_int(n_act_orb,2) double precision :: hij,phase !matrix_2h1p = 0.d0 @@ -34,8 +34,8 @@ subroutine give_2h1p_contrib(matrix_2h1p) ! take all the integral you will need for i,j,r fixed do a = 1, n_act_orb aorb = list_act(a) - active_int(a,1) = get_mo_bielec_integral_schwartz(iorb,jorb,rorb,aorb,mo_integrals_map) ! direct - active_int(a,2) = get_mo_bielec_integral_schwartz(iorb,jorb,aorb,rorb,mo_integrals_map) ! exchange + active_int(a,1) = get_mo_bielec_integral(iorb,jorb,rorb,aorb,mo_integrals_map) ! direct + active_int(a,2) = get_mo_bielec_integral(iorb,jorb,aorb,rorb,mo_integrals_map) ! exchange enddo integer :: degree(N_det) @@ -209,7 +209,7 @@ subroutine give_1h2p_contrib(matrix_1h2p) integer(bit_kind) :: det_tmp(N_int,2) integer :: exc(0:2,2,2) integer :: accu_elec - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral double precision :: active_int(n_act_orb,2) double precision :: hij,phase !matrix_1h2p = 0.d0 @@ -228,8 +228,8 @@ subroutine give_1h2p_contrib(matrix_1h2p) ! take all the integral you will need for i,j,r fixed do a = 1, n_act_orb aorb = list_act(a) - active_int(a,1) = get_mo_bielec_integral_schwartz(iorb,aorb,rorb,vorb,mo_integrals_map) ! direct - active_int(a,2) = get_mo_bielec_integral_schwartz(iorb,aorb,vorb,rorb,mo_integrals_map) ! exchange + active_int(a,1) = get_mo_bielec_integral(iorb,aorb,rorb,vorb,mo_integrals_map) ! direct + active_int(a,2) = get_mo_bielec_integral(iorb,aorb,vorb,rorb,mo_integrals_map) ! exchange enddo integer :: degree(N_det) @@ -406,7 +406,7 @@ subroutine give_1h1p_contrib(matrix_1h1p) integer(bit_kind) :: det_tmp(N_int,2) integer :: exc(0:2,2,2) integer :: accu_elec - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral double precision :: active_int(n_act_orb,2) double precision :: hij,phase integer :: degree(N_det) @@ -474,10 +474,10 @@ subroutine give_1h1p_contrib(matrix_1h1p) endif call get_double_excitation(psi_det(1,1,idx(jdet)),det_tmp,exc,phase,N_int) if(ispin == jspin )then - hij = -get_mo_bielec_integral_schwartz(iorb,aorb,rorb,borb,mo_integrals_map) & - + get_mo_bielec_integral_schwartz(iorb,aorb,borb,rorb,mo_integrals_map) + hij = -get_mo_bielec_integral(iorb,aorb,rorb,borb,mo_integrals_map) & + + get_mo_bielec_integral(iorb,aorb,borb,rorb,mo_integrals_map) else - hij = get_mo_bielec_integral_schwartz(iorb,borb,rorb,aorb,mo_integrals_map) + hij = get_mo_bielec_integral(iorb,borb,rorb,aorb,mo_integrals_map) endif hij = hij * phase double precision :: hij_test @@ -530,7 +530,7 @@ subroutine give_1h1p_sec_order_singles_contrib(matrix_1h1p) double precision :: hij_det_pert(n_inact_orb,n_virt_orb,2,N_states) integer :: exc(0:2,2,2) integer :: accu_elec - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral double precision :: active_int(n_act_orb,2) double precision :: hij,phase integer :: degree(N_det) @@ -690,7 +690,7 @@ subroutine give_1p_sec_order_singles_contrib(matrix_1p) double precision :: hij_det_pert(n_act_orb,n_virt_orb,2) integer :: exc(0:2,2,2) integer :: accu_elec - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral double precision :: hij,phase integer :: degree(N_det) integer :: idx(0:N_det) @@ -832,7 +832,7 @@ subroutine give_1h1p_only_doubles_spin_cross(matrix_1h1p) integer(bit_kind) :: det_tmp(N_int,2) integer :: exc(0:2,2,2) integer :: accu_elec - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral double precision :: active_int(n_act_orb,2) double precision :: hij,phase integer :: degree(N_det) diff --git a/plugins/MRPT_Utils/new_way_second_order_coef.irp.f b/plugins/MRPT_Utils/new_way_second_order_coef.irp.f index 5c4b562f..4c12dbe1 100644 --- a/plugins/MRPT_Utils/new_way_second_order_coef.irp.f +++ b/plugins/MRPT_Utils/new_way_second_order_coef.irp.f @@ -16,7 +16,7 @@ subroutine give_2h1p_contrib_sec_order(matrix_2h1p) integer(bit_kind) :: det_tmp_j(N_int,2) integer :: exc(0:2,2,2) integer :: accu_elec - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral double precision :: active_int(n_act_orb,2) double precision :: hij,phase integer :: index_orb_act_mono(N_det,6) @@ -36,8 +36,8 @@ subroutine give_2h1p_contrib_sec_order(matrix_2h1p) ! take all the integral you will need for i,j,r fixed do a = 1, n_act_orb aorb = list_act(a) - active_int(a,1) = get_mo_bielec_integral_schwartz(iorb,jorb,rorb,aorb,mo_integrals_map) ! direct - active_int(a,2) = get_mo_bielec_integral_schwartz(iorb,jorb,aorb,rorb,mo_integrals_map) ! exchange + active_int(a,1) = get_mo_bielec_integral(iorb,jorb,rorb,aorb,mo_integrals_map) ! direct + active_int(a,2) = get_mo_bielec_integral(iorb,jorb,aorb,rorb,mo_integrals_map) ! exchange perturb_dets_phase(a,1,1) = -1000.d0 perturb_dets_phase(a,1,2) = -1000.d0 perturb_dets_phase(a,2,2) = -1000.d0 @@ -375,7 +375,7 @@ subroutine give_1h2p_contrib_sec_order(matrix_1h2p) integer(bit_kind) :: det_tmp_j(N_int,2) integer :: exc(0:2,2,2) integer :: accu_elec - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral double precision :: active_int(n_act_orb,2) double precision :: hij,phase double precision :: accu_contrib @@ -410,8 +410,8 @@ subroutine give_1h2p_contrib_sec_order(matrix_1h2p) ! take all the integral you will need for i,j,r fixed do a = 1, n_act_orb aorb = list_act(a) - active_int(a,1) = get_mo_bielec_integral_schwartz(iorb,aorb,rorb,vorb,mo_integrals_map) ! direct - active_int(a,2) = get_mo_bielec_integral_schwartz(iorb,aorb,vorb,rorb,mo_integrals_map) ! exchange + active_int(a,1) = get_mo_bielec_integral(iorb,aorb,rorb,vorb,mo_integrals_map) ! direct + active_int(a,2) = get_mo_bielec_integral(iorb,aorb,vorb,rorb,mo_integrals_map) ! exchange perturb_dets_phase(a,1,1) = -1000.d0 perturb_dets_phase(a,1,2) = -1000.d0 perturb_dets_phase(a,2,2) = -1000.d0 diff --git a/plugins/MRPT_Utils/second_order_new.irp.f b/plugins/MRPT_Utils/second_order_new.irp.f index bcd08bf5..46de6601 100644 --- a/plugins/MRPT_Utils/second_order_new.irp.f +++ b/plugins/MRPT_Utils/second_order_new.irp.f @@ -18,7 +18,7 @@ subroutine give_1h2p_new(matrix_1h2p) integer(bit_kind) :: det_tmp_j(N_int,2) integer :: exc(0:2,2,2) integer :: accu_elec - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral double precision :: active_int(n_act_orb,2) double precision :: hij,phase double precision :: accu_contrib(N_states) @@ -63,8 +63,8 @@ subroutine give_1h2p_new(matrix_1h2p) ! take all the integral you will need for i,j,r fixed do a = 1, n_act_orb aorb = list_act(a) - active_int(a,1) = get_mo_bielec_integral_schwartz(iorb,aorb,rorb,vorb,mo_integrals_map) ! direct - active_int(a,2) = get_mo_bielec_integral_schwartz(iorb,aorb,vorb,rorb,mo_integrals_map) ! exchange + active_int(a,1) = get_mo_bielec_integral(iorb,aorb,rorb,vorb,mo_integrals_map) ! direct + active_int(a,2) = get_mo_bielec_integral(iorb,aorb,vorb,rorb,mo_integrals_map) ! exchange perturb_dets_phase(a,1,1) = -1000.d0 perturb_dets_phase(a,1,2) = -1000.d0 perturb_dets_phase(a,2,2) = -1000.d0 @@ -495,7 +495,7 @@ subroutine give_2h1p_new(matrix_2h1p) integer(bit_kind) :: det_tmp(N_int,2) integer :: exc(0:2,2,2) integer :: accu_elec - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral double precision :: active_int(n_act_orb,2) double precision :: hij,phase integer :: i_hole,i_part @@ -531,8 +531,8 @@ subroutine give_2h1p_new(matrix_2h1p) ! take all the integral you will need for i,j,r fixed do a = 1, n_act_orb aorb = list_act(a) - active_int(a,1) = get_mo_bielec_integral_schwartz(iorb,jorb,rorb,aorb,mo_integrals_map) ! direct - active_int(a,2) = get_mo_bielec_integral_schwartz(iorb,jorb,aorb,rorb,mo_integrals_map) ! exchange + active_int(a,1) = get_mo_bielec_integral(iorb,jorb,rorb,aorb,mo_integrals_map) ! direct + active_int(a,2) = get_mo_bielec_integral(iorb,jorb,aorb,rorb,mo_integrals_map) ! exchange perturb_dets_phase(a,1,1) = -1000.d0 perturb_dets_phase(a,1,2) = -1000.d0 perturb_dets_phase(a,2,2) = -1000.d0 diff --git a/plugins/MRPT_Utils/second_order_new_2p.irp.f b/plugins/MRPT_Utils/second_order_new_2p.irp.f index 2e94527c..11ae18da 100644 --- a/plugins/MRPT_Utils/second_order_new_2p.irp.f +++ b/plugins/MRPT_Utils/second_order_new_2p.irp.f @@ -17,7 +17,7 @@ subroutine give_2p_new(matrix_2p) integer(bit_kind) :: det_tmp_j(N_int,2) integer :: exc(0:2,2,2) integer :: accu_elec - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral double precision :: active_int(n_act_orb,n_act_orb,2) double precision :: hij,phase double precision :: accu_contrib(N_states) @@ -62,8 +62,8 @@ subroutine give_2p_new(matrix_2p) aorb = list_act(a) do b = 1, n_act_orb borb = list_act(b) - active_int(a,b,1) = get_mo_bielec_integral_schwartz(aorb,borb,rorb,vorb,mo_integrals_map) ! direct ( a--> r | b--> v ) - active_int(a,b,2) = get_mo_bielec_integral_schwartz(aorb,borb,vorb,rorb,mo_integrals_map) ! exchange ( b--> r | a--> v ) + active_int(a,b,1) = get_mo_bielec_integral(aorb,borb,rorb,vorb,mo_integrals_map) ! direct ( a--> r | b--> v ) + active_int(a,b,2) = get_mo_bielec_integral(aorb,borb,vorb,rorb,mo_integrals_map) ! exchange ( b--> r | a--> v ) perturb_dets_phase(a,b,1,1) = -1000.d0 perturb_dets_phase(a,b,1,2) = -1000.d0 perturb_dets_phase(a,b,2,2) = -1000.d0 diff --git a/plugins/Perturbation/EZFIO.cfg b/plugins/Perturbation/EZFIO.cfg index c5d6379f..9023accf 100644 --- a/plugins/Perturbation/EZFIO.cfg +++ b/plugins/Perturbation/EZFIO.cfg @@ -17,3 +17,4 @@ doc: The selection process stops when the energy ratio variational/(variational+ is equal to var_pt2_ratio interface: ezfio,provider,ocaml default: 0.75 + diff --git a/plugins/Properties/EZFIO.cfg b/plugins/Properties/EZFIO.cfg index 40ccd8b9..2a5ae803 100644 --- a/plugins/Properties/EZFIO.cfg +++ b/plugins/Properties/EZFIO.cfg @@ -9,3 +9,4 @@ type: double precision doc: threshold for the values of the alpha/beta two body dm evaluation interface: ezfio,provider,ocaml default: 0.000001 + diff --git a/plugins/Properties/hyperfine_constants.irp.f b/plugins/Properties/hyperfine_constants.irp.f index 63ad545b..6fa39278 100644 --- a/plugins/Properties/hyperfine_constants.irp.f +++ b/plugins/Properties/hyperfine_constants.irp.f @@ -151,7 +151,7 @@ subroutine print_hcc integer :: i,j print*,'Z AU GAUSS MHZ cm^-1' do i = 1, nucl_num - write(*,'(I2,X,F3.1,X,4(F16.6,X))')i,nucl_charge(i),spin_density_at_nucleous(i),iso_hcc_gauss(i),iso_hcc_mhz(i),iso_hcc_cm_1(i) + write(*,'(I2,X,F4.1,X,4(F16.6,X))')i,nucl_charge(i),spin_density_at_nucleous(i),iso_hcc_gauss(i),iso_hcc_mhz(i),iso_hcc_cm_1(i) enddo end diff --git a/plugins/Selectors_full/selectors.irp.f b/plugins/Selectors_full/selectors.irp.f index 71c3550e..ce5e8367 100644 --- a/plugins/Selectors_full/selectors.irp.f +++ b/plugins/Selectors_full/selectors.irp.f @@ -20,8 +20,7 @@ BEGIN_PROVIDER [ integer, N_det_selectors] norm = norm + psi_average_norm_contrib_sorted(i) if (norm > threshold_selectors) then -! N_det_selectors = i-1 - N_det_selectors = i + N_det_selectors = i-1 exit endif enddo diff --git a/src/Bitmask/bitmask_cas_routines.irp.f b/src/Bitmask/bitmask_cas_routines.irp.f index 961be5df..87a02d10 100644 --- a/src/Bitmask/bitmask_cas_routines.irp.f +++ b/src/Bitmask/bitmask_cas_routines.irp.f @@ -1,7 +1,8 @@ use bitmasks integer function number_of_holes(key_in) -use bitmasks - ! function that returns the number of holes in the inact space + BEGIN_DOC + ! Function that returns the number of holes in the inact space + END_DOC implicit none integer(bit_kind), intent(in) :: key_in(N_int,2) integer :: i @@ -104,8 +105,9 @@ end integer function number_of_particles(key_in) -use bitmasks + BEGIN_DOC ! function that returns the number of particles in the virtual space + END_DOC implicit none integer(bit_kind), intent(in) :: key_in(N_int,2) integer :: i @@ -208,12 +210,13 @@ use bitmasks end logical function is_a_two_holes_two_particles(key_in) -use bitmasks + BEGIN_DOC ! logical function that returns True if the determinant 'key_in' ! belongs to the 2h-2p excitation class of the DDCI space ! this is calculated using the CAS_bitmask that defines the active ! orbital space, the inact_bitmasl that defines the inactive oribital space ! and the virt_bitmask that defines the virtual orbital space + END_DOC implicit none integer(bit_kind), intent(in) :: key_in(N_int,2) integer :: i,i_diff @@ -403,8 +406,9 @@ use bitmasks integer function number_of_holes_verbose(key_in) -use bitmasks + BEGIN_DOC ! function that returns the number of holes in the inact space + END_DOC implicit none integer(bit_kind), intent(in) :: key_in(N_int,2) integer :: i @@ -432,7 +436,9 @@ end integer function number_of_particles_verbose(key_in) + BEGIN_DOC ! function that returns the number of particles in the inact space + END_DOC implicit none integer(bit_kind), intent(in) :: key_in(N_int,2) integer :: i @@ -458,7 +464,6 @@ integer function number_of_particles_verbose(key_in) end logical function is_a_1h1p(key_in) -use bitmasks implicit none integer(bit_kind), intent(in) :: key_in(N_int,2) integer :: number_of_particles, number_of_holes @@ -470,7 +475,6 @@ use bitmasks end logical function is_a_1h2p(key_in) -use bitmasks implicit none integer(bit_kind), intent(in) :: key_in(N_int,2) integer :: number_of_particles, number_of_holes @@ -482,7 +486,6 @@ use bitmasks end logical function is_a_2h1p(key_in) -use bitmasks implicit none integer(bit_kind), intent(in) :: key_in(N_int,2) integer :: number_of_particles, number_of_holes @@ -494,7 +497,6 @@ use bitmasks end logical function is_a_1h(key_in) -use bitmasks implicit none integer(bit_kind), intent(in) :: key_in(N_int,2) integer :: number_of_particles, number_of_holes @@ -506,7 +508,6 @@ use bitmasks end logical function is_a_1p(key_in) -use bitmasks implicit none integer(bit_kind), intent(in) :: key_in(N_int,2) integer :: number_of_particles, number_of_holes @@ -518,7 +519,6 @@ use bitmasks end logical function is_a_2p(key_in) -use bitmasks implicit none integer(bit_kind), intent(in) :: key_in(N_int,2) integer :: number_of_particles, number_of_holes @@ -530,7 +530,6 @@ use bitmasks end logical function is_a_2h(key_in) -use bitmasks implicit none integer(bit_kind), intent(in) :: key_in(N_int,2) integer :: number_of_particles, number_of_holes @@ -542,7 +541,6 @@ use bitmasks end logical function is_i_in_virtual(i) -use bitmasks implicit none integer,intent(in) :: i integer(bit_kind) :: key(N_int) diff --git a/src/Determinants/EZFIO.cfg b/src/Determinants/EZFIO.cfg index 324005aa..afb2644e 100644 --- a/src/Determinants/EZFIO.cfg +++ b/src/Determinants/EZFIO.cfg @@ -127,13 +127,13 @@ default: 0. [store_full_H_mat] type: logical -doc: If True, the Davidson diagonalization is performed by storring the full H matrix up to n_det_max_stored. Be carefull, it can cost a lot of memory but can also save a lot of CPU time +doc: If True, the Davidson diagonalization is performed by storing the full H matrix up to n_det_max_stored. Be careful, it can cost a lot of memory but can also save a lot of CPU time interface: ezfio,provider,ocaml default: False [n_det_max_stored] type: Det_number_max -doc: Maximum number of determinants for which the full H matrix is stored. Be carefull, the memory requested scales as 10*n_det_max_stored**2. For instance, 90000 determinants represent a matrix of size 60 Gb. +doc: Maximum number of determinants for which the full H matrix is stored. Be careful, the memory requested scales as 10*n_det_max_stored**2. For instance, 90000 determinants represent a matrix of size 60 Gb. interface: ezfio,provider,ocaml default: 90000 diff --git a/src/Determinants/NEEDED_CHILDREN_MODULES b/src/Determinants/NEEDED_CHILDREN_MODULES index 566762ba..8711010f 100644 --- a/src/Determinants/NEEDED_CHILDREN_MODULES +++ b/src/Determinants/NEEDED_CHILDREN_MODULES @@ -1 +1 @@ -Integrals_Monoelec Integrals_Bielec Hartree_Fock +Integrals_Monoelec Integrals_Bielec diff --git a/src/Determinants/diagonalize_CI.irp.f b/src/Determinants/diagonalize_CI.irp.f index 49714082..6f94eedb 100644 --- a/src/Determinants/diagonalize_CI.irp.f +++ b/src/Determinants/diagonalize_CI.irp.f @@ -77,14 +77,10 @@ END_PROVIDER if (diag_algorithm == "Davidson") then - print*, '------------- In Davidson ' call davidson_diag(psi_det,CI_eigenvectors,CI_electronic_energy, & size(CI_eigenvectors,1),N_det,N_states_diag,N_int,output_determinants) - print*, '------------- Out Davidson ' do j=1,N_states_diag - print*, '------------- In S^2' call get_s2_u0(psi_det,CI_eigenvectors(1,j),N_det,size(CI_eigenvectors,1),CI_eigenvectors_s2(j)) - print*, '------------- Out S^2' enddo @@ -103,7 +99,6 @@ END_PROVIDER do j=1,N_det call get_s2_u0(psi_det,eigenvectors(1,j),N_det,size(eigenvectors,1),s2) s2_eigvalues(j) = s2 - print*, 's2 in lapack',s2 print*, eigenvalues(j) + nuclear_repulsion ! Select at least n_states states with S^2 values closed to "expected_s2" if(dabs(s2-expected_s2).le.0.3d0)then @@ -219,12 +214,6 @@ END_PROVIDER do i = 1, N_det CI_eigenvectors(i,j) = psi_coef(i,index_good_state_array(iorder(j))) enddo -! call u0_H_u_0(e_0,CI_eigenvectors(1,j),n_det,psi_det,N_int) -! print*,'e = ',CI_electronic_energy(j) -! print*,' = ',e_0 -! call get_s2_u0(psi_det,CI_eigenvectors(1,j),N_det,size(CI_eigenvectors,1),s2) -! print*,'s^2 = ',CI_eigenvectors_s2(j) -! print*,'= ',s2 enddo deallocate(e_array,iorder) @@ -269,7 +258,6 @@ END_PROVIDER endif deallocate(s2_eigvalues) endif - print*, 'out provider' END_PROVIDER diff --git a/src/Determinants/slater_rules.irp.f b/src/Determinants/slater_rules.irp.f index 5cbed15e..f98947a2 100644 --- a/src/Determinants/slater_rules.irp.f +++ b/src/Determinants/slater_rules.irp.f @@ -443,7 +443,7 @@ subroutine i_H_j_new(key_i,key_j,Nint,hij) integer :: exc(0:2,2,2) integer :: degree - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral integer :: m,n,p,q integer :: i,j,k integer :: occ(Nint*bit_kind_size,2) @@ -469,31 +469,31 @@ subroutine i_H_j_new(key_i,key_j,Nint,hij) call get_double_excitation(key_i,key_j,exc,phase,Nint) if (exc(0,1,1) == 1) then ! Mono alpha, mono beta - hij = phase*get_mo_bielec_integral_schwartz( & + hij = phase*get_mo_bielec_integral( & exc(1,1,1), & exc(1,1,2), & exc(1,2,1), & exc(1,2,2) ,mo_integrals_map) else if (exc(0,1,1) == 2) then ! Double alpha - hij = phase*(get_mo_bielec_integral_schwartz( & + hij = phase*(get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(1,2,1), & exc(2,2,1) ,mo_integrals_map) - & - get_mo_bielec_integral_schwartz( & + get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(2,2,1), & exc(1,2,1) ,mo_integrals_map) ) else if (exc(0,1,2) == 2) then ! Double beta - hij = phase*(get_mo_bielec_integral_schwartz( & + hij = phase*(get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(1,2,2), & exc(2,2,2) ,mo_integrals_map) - & - get_mo_bielec_integral_schwartz( & + get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(2,2,2), & @@ -512,15 +512,15 @@ subroutine i_H_j_new(key_i,key_j,Nint,hij) do k = 1, elec_alpha_num i = occ(k,1) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) - miip(i) = get_mo_bielec_integral_schwartz(m,i,i,p,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) + miip(i) = get_mo_bielec_integral(m,i,i,p,mo_integrals_map) has_mipi(i) = .True. endif enddo do k = 1, elec_beta_num i = occ(k,2) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) has_mipi(i) = .True. endif enddo @@ -541,15 +541,15 @@ subroutine i_H_j_new(key_i,key_j,Nint,hij) do k = 1, elec_beta_num i = occ(k,2) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) - miip(i) = get_mo_bielec_integral_schwartz(m,i,i,p,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) + miip(i) = get_mo_bielec_integral(m,i,i,p,mo_integrals_map) has_mipi(i) = .True. endif enddo do k = 1, elec_alpha_num i = occ(k,1) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) has_mipi(i) = .True. endif enddo @@ -583,7 +583,7 @@ subroutine i_H_j(key_i,key_j,Nint,hij) integer :: exc(0:2,2,2) integer :: degree - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral integer :: m,n,p,q integer :: i,j,k integer :: occ(Nint*bit_kind_size,2) @@ -614,7 +614,7 @@ subroutine i_H_j(key_i,key_j,Nint,hij) else if (exc(1,2,1) ==exc(1,1,2))then hij = phase * big_array_exchange_integrals(exc(1,2,1),exc(1,1,1),exc(1,2,2)) else - hij = phase*get_mo_bielec_integral_schwartz( & + hij = phase*get_mo_bielec_integral( & exc(1,1,1), & exc(1,1,2), & exc(1,2,1), & @@ -622,24 +622,24 @@ subroutine i_H_j(key_i,key_j,Nint,hij) endif else if (exc(0,1,1) == 2) then ! Double alpha - hij = phase*(get_mo_bielec_integral_schwartz( & + hij = phase*(get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(1,2,1), & exc(2,2,1) ,mo_integrals_map) - & - get_mo_bielec_integral_schwartz( & + get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(2,2,1), & exc(1,2,1) ,mo_integrals_map) ) else if (exc(0,1,2) == 2) then ! Double beta - hij = phase*(get_mo_bielec_integral_schwartz( & + hij = phase*(get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(1,2,2), & exc(2,2,2) ,mo_integrals_map) - & - get_mo_bielec_integral_schwartz( & + get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(2,2,2), & @@ -658,15 +658,15 @@ subroutine i_H_j(key_i,key_j,Nint,hij) ! do k = 1, elec_alpha_num ! i = occ(k,1) ! if (.not.has_mipi(i)) then -! mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) -! miip(i) = get_mo_bielec_integral_schwartz(m,i,i,p,mo_integrals_map) +! mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) +! miip(i) = get_mo_bielec_integral(m,i,i,p,mo_integrals_map) ! has_mipi(i) = .True. ! endif ! enddo ! do k = 1, elec_beta_num ! i = occ(k,2) ! if (.not.has_mipi(i)) then -! mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) +! mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) ! has_mipi(i) = .True. ! endif ! enddo @@ -687,15 +687,15 @@ subroutine i_H_j(key_i,key_j,Nint,hij) ! do k = 1, elec_beta_num ! i = occ(k,2) ! if (.not.has_mipi(i)) then -! mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) -! miip(i) = get_mo_bielec_integral_schwartz(m,i,i,p,mo_integrals_map) +! mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) +! miip(i) = get_mo_bielec_integral(m,i,i,p,mo_integrals_map) ! has_mipi(i) = .True. ! endif ! enddo ! do k = 1, elec_alpha_num ! i = occ(k,1) ! if (.not.has_mipi(i)) then -! mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) +! mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) ! has_mipi(i) = .True. ! endif ! enddo @@ -731,7 +731,7 @@ subroutine i_H_j_phase_out(key_i,key_j,Nint,hij,phase,exc,degree) integer,intent(out) :: exc(0:2,2,2) integer,intent(out) :: degree - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral integer :: m,n,p,q integer :: i,j,k integer :: occ(Nint*bit_kind_size,2) @@ -756,31 +756,31 @@ subroutine i_H_j_phase_out(key_i,key_j,Nint,hij,phase,exc,degree) call get_double_excitation(key_i,key_j,exc,phase,Nint) if (exc(0,1,1) == 1) then ! Mono alpha, mono beta - hij = phase*get_mo_bielec_integral_schwartz( & + hij = phase*get_mo_bielec_integral( & exc(1,1,1), & exc(1,1,2), & exc(1,2,1), & exc(1,2,2) ,mo_integrals_map) else if (exc(0,1,1) == 2) then ! Double alpha - hij = phase*(get_mo_bielec_integral_schwartz( & + hij = phase*(get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(1,2,1), & exc(2,2,1) ,mo_integrals_map) - & - get_mo_bielec_integral_schwartz( & + get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(2,2,1), & exc(1,2,1) ,mo_integrals_map) ) else if (exc(0,1,2) == 2) then ! Double beta - hij = phase*(get_mo_bielec_integral_schwartz( & + hij = phase*(get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(1,2,2), & exc(2,2,2) ,mo_integrals_map) - & - get_mo_bielec_integral_schwartz( & + get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(2,2,2), & @@ -798,15 +798,15 @@ subroutine i_H_j_phase_out(key_i,key_j,Nint,hij,phase,exc,degree) do k = 1, elec_alpha_num i = occ(k,1) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) - miip(i) = get_mo_bielec_integral_schwartz(m,i,i,p,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) + miip(i) = get_mo_bielec_integral(m,i,i,p,mo_integrals_map) has_mipi(i) = .True. endif enddo do k = 1, elec_beta_num i = occ(k,2) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) has_mipi(i) = .True. endif enddo @@ -825,15 +825,15 @@ subroutine i_H_j_phase_out(key_i,key_j,Nint,hij,phase,exc,degree) do k = 1, elec_beta_num i = occ(k,2) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) - miip(i) = get_mo_bielec_integral_schwartz(m,i,i,p,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) + miip(i) = get_mo_bielec_integral(m,i,i,p,mo_integrals_map) has_mipi(i) = .True. endif enddo do k = 1, elec_alpha_num i = occ(k,1) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) has_mipi(i) = .True. endif enddo @@ -867,7 +867,7 @@ subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble) integer :: exc(0:2,2,2) integer :: degree - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral integer :: m,n,p,q integer :: i,j,k integer :: occ(Nint*bit_kind_size,2) @@ -894,7 +894,7 @@ subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble) call get_double_excitation(key_i,key_j,exc,phase,Nint) if (exc(0,1,1) == 1) then ! Mono alpha, mono beta - hij = phase*get_mo_bielec_integral_schwartz( & + hij = phase*get_mo_bielec_integral( & exc(1,1,1), & exc(1,1,2), & exc(1,2,1), & @@ -904,22 +904,22 @@ subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble) else if (exc(0,1,1) == 2) then ! Double alpha print*,'phase hij = ',phase - hij = phase*(get_mo_bielec_integral_schwartz( & + hij = phase*(get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(1,2,1), & exc(2,2,1) ,mo_integrals_map) - & - get_mo_bielec_integral_schwartz( & + get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(2,2,1), & exc(1,2,1) ,mo_integrals_map) ) - print*,get_mo_bielec_integral_schwartz( & + print*,get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(1,2,1), & exc(2,2,1) ,mo_integrals_map) - print*,get_mo_bielec_integral_schwartz( & + print*,get_mo_bielec_integral( & exc(1,1,1), & exc(2,1,1), & exc(2,2,1), & @@ -928,23 +928,23 @@ subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble) else if (exc(0,1,2) == 2) then ! Double beta print*,'phase hij = ',phase - print*, get_mo_bielec_integral_schwartz( & + print*, get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(1,2,2), & exc(2,2,2) ,mo_integrals_map ) - print*, get_mo_bielec_integral_schwartz( & + print*, get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(2,2,2), & exc(1,2,2) ,mo_integrals_map) - hij = phase*(get_mo_bielec_integral_schwartz( & + hij = phase*(get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(1,2,2), & exc(2,2,2) ,mo_integrals_map) - & - get_mo_bielec_integral_schwartz( & + get_mo_bielec_integral( & exc(1,1,2), & exc(2,1,2), & exc(2,2,2), & @@ -962,15 +962,15 @@ subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble) do k = 1, elec_alpha_num i = occ(k,1) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) - miip(i) = get_mo_bielec_integral_schwartz(m,i,i,p,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) + miip(i) = get_mo_bielec_integral(m,i,i,p,mo_integrals_map) has_mipi(i) = .True. endif enddo do k = 1, elec_beta_num i = occ(k,2) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) has_mipi(i) = .True. endif enddo @@ -989,15 +989,15 @@ subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble) do k = 1, elec_beta_num i = occ(k,2) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) - miip(i) = get_mo_bielec_integral_schwartz(m,i,i,p,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) + miip(i) = get_mo_bielec_integral(m,i,i,p,mo_integrals_map) has_mipi(i) = .True. endif enddo do k = 1, elec_alpha_num i = occ(k,1) if (.not.has_mipi(i)) then - mipi(i) = get_mo_bielec_integral_schwartz(m,i,p,i,mo_integrals_map) + mipi(i) = get_mo_bielec_integral(m,i,p,i,mo_integrals_map) has_mipi(i) = .True. endif enddo diff --git a/src/Determinants/test_3d.irp.f b/src/Determinants/test_3d.irp.f deleted file mode 100644 index a5d09cd3..00000000 --- a/src/Determinants/test_3d.irp.f +++ /dev/null @@ -1,15 +0,0 @@ -program test_3d - implicit none - integer :: i,npt - double precision :: dx,domain,x_min,x,step_function_becke -!domain = 5.d0 -!npt = 100 -!dx = domain/dble(npt) -!x_min = -0.5d0 * domain -!x = x_min -!do i = 1, npt -! write(33,*)x,step_function_becke(x) -! x += dx -!enddo - -end diff --git a/src/Determinants/test_two_body.irp.f b/src/Determinants/test_two_body.irp.f deleted file mode 100644 index 54c43c09..00000000 --- a/src/Determinants/test_two_body.irp.f +++ /dev/null @@ -1,18 +0,0 @@ -program test - implicit none - read_wf = .True. - touch read_wf - call routine -end - -subroutine routine - implicit none - integer :: i,j,k,l - do i = 1, n_act_orb - do j = 1, n_act_orb - do k = 1, n_act_orb - - enddo - enddo - enddo -end diff --git a/src/Determinants/utils.irp.f b/src/Determinants/utils.irp.f index cc191970..dbd5a7ef 100644 --- a/src/Determinants/utils.irp.f +++ b/src/Determinants/utils.irp.f @@ -8,13 +8,9 @@ BEGIN_PROVIDER [ double precision, H_matrix_all_dets,(N_det,N_det) ] double precision :: hij integer :: degree(N_det),idx(0:N_det) call i_H_j(psi_det(1,1,1),psi_det(1,1,1),N_int,hij) - !$OMP PARALLEL DO SCHEDULE(GUIDED) PRIVATE(i,j,hij,degree,idx,k) & + !$OMP PARALLEL DO SCHEDULE(GUIDED) DEFAULT(NONE) PRIVATE(i,j,hij,degree,idx,k) & !$OMP SHARED (N_det, psi_det, N_int,H_matrix_all_dets) do i =1,N_det -! call get_excitation_degree_vector(psi_det,psi_det(1,1,i),degree,N_int,N_det,idx) -! do k =1, idx(0) -! j = idx(k) -! if(j.lt.i)cycle do j = i, N_det call i_H_j(psi_det(1,1,i),psi_det(1,1,j),N_int,hij) H_matrix_all_dets(i,j) = hij @@ -25,32 +21,3 @@ BEGIN_PROVIDER [ double precision, H_matrix_all_dets,(N_det,N_det) ] END_PROVIDER -subroutine provide_big_matrix_stored_with_current_dets(sze,dets_in,big_matrix_stored) - use bitmasks - integer, intent(in) :: sze - integer(bit_kind), intent(in) :: dets_in(N_int,2,sze) - double precision, intent(out) :: big_matrix_stored(sze,sze) - integer :: i,j,k - double precision :: hij - integer :: degree(N_det),idx(0:N_det) - call i_H_j(dets_in(1,1,1),dets_in(1,1,1),N_int,hij) - print*, 'providing big_matrix_stored' - print*, n_det_max_stored - !$OMP PARALLEL DO SCHEDULE(GUIDED) PRIVATE(i,j,hij,degree,idx,k) & - !$OMP SHARED (sze, dets_in, N_int,big_matrix_stored) - do i =1,sze -! call get_excitation_degree_vector(dets_in,dets_in(1,1,i),degree,N_int,sze,idx) -! do k =1, idx(0) -! j = idx(k) - do j = i, sze - if(j.lt.i)cycle - call i_H_j(dets_in(1,1,i),dets_in(1,1,j),N_int,hij) - big_matrix_stored(i,j) = hij - big_matrix_stored(j,i) = hij - enddo - enddo - !$OMP END PARALLEL DO - print*, 'big_matrix_stored provided !!' - - -end diff --git a/src/Integrals_Bielec/integrals_3_index.irp.f b/src/Integrals_Bielec/integrals_3_index.irp.f index b9ee29b9..41037b34 100644 --- a/src/Integrals_Bielec/integrals_3_index.irp.f +++ b/src/Integrals_Bielec/integrals_3_index.irp.f @@ -2,17 +2,17 @@ &BEGIN_PROVIDER [double precision, big_array_exchange_integrals,(mo_tot_num_align,mo_tot_num, mo_tot_num)] implicit none integer :: i,j,k,l - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral double precision :: integral do k = 1, mo_tot_num do i = 1, mo_tot_num do j = 1, mo_tot_num l = j - integral = get_mo_bielec_integral_schwartz(i,j,k,l,mo_integrals_map) + integral = get_mo_bielec_integral(i,j,k,l,mo_integrals_map) big_array_coulomb_integrals(j,i,k) = integral l = j - integral = get_mo_bielec_integral_schwartz(i,j,l,k,mo_integrals_map) + integral = get_mo_bielec_integral(i,j,l,k,mo_integrals_map) big_array_exchange_integrals(j,i,k) = integral enddo enddo diff --git a/src/Integrals_Bielec/map_integrals.irp.f b/src/Integrals_Bielec/map_integrals.irp.f index 305abee3..65561a57 100644 --- a/src/Integrals_Bielec/map_integrals.irp.f +++ b/src/Integrals_Bielec/map_integrals.irp.f @@ -294,28 +294,6 @@ double precision function get_mo_bielec_integral(i,j,k,l,map) get_mo_bielec_integral = dble(tmp) end -double precision function get_mo_bielec_integral_schwartz(i,j,k,l,map) - use map_module - implicit none - BEGIN_DOC - ! Returns one integral in the MO basis - END_DOC - integer, intent(in) :: i,j,k,l - integer(key_kind) :: idx - type(map_type), intent(inout) :: map - real(integral_kind) :: tmp - PROVIDE mo_bielec_integrals_in_map - if (mo_bielec_integral_schwartz(i,k)*mo_bielec_integral_schwartz(j,l) > mo_integrals_threshold) then - !DIR$ FORCEINLINE - call bielec_integrals_index(i,j,k,l,idx) - !DIR$ FORCEINLINE - call map_get(map,idx,tmp) - else - tmp = 0.d0 - endif - get_mo_bielec_integral_schwartz = dble(tmp) -end - double precision function mo_bielec_integral(i,j,k,l) implicit none @@ -323,9 +301,9 @@ double precision function mo_bielec_integral(i,j,k,l) ! Returns one integral in the MO basis END_DOC integer, intent(in) :: i,j,k,l - double precision :: get_mo_bielec_integral_schwartz + double precision :: get_mo_bielec_integral PROVIDE mo_bielec_integrals_in_map - mo_bielec_integral = get_mo_bielec_integral_schwartz(i,j,k,l,mo_integrals_map) + mo_bielec_integral = get_mo_bielec_integral(i,j,k,l,mo_integrals_map) return end diff --git a/src/Integrals_Bielec/mo_bi_integrals.irp.f b/src/Integrals_Bielec/mo_bi_integrals.irp.f index bf23ad1f..af7f21d2 100644 --- a/src/Integrals_Bielec/mo_bi_integrals.irp.f +++ b/src/Integrals_Bielec/mo_bi_integrals.irp.f @@ -1447,22 +1447,6 @@ END_PROVIDER END_PROVIDER -BEGIN_PROVIDER [ double precision, mo_bielec_integral_schwartz,(mo_tot_num,mo_tot_num) ] - implicit none - BEGIN_DOC - ! Needed to compute Schwartz inequalities - END_DOC - - integer :: i,k - - do i=1,mo_tot_num - do k=1,mo_tot_num - mo_bielec_integral_schwartz(k,i) = 1.d10 - enddo - enddo - -END_PROVIDER - subroutine clear_mo_map implicit none @@ -1470,7 +1454,7 @@ subroutine clear_mo_map ! Frees the memory of the MO map END_DOC call map_deinit(mo_integrals_map) - FREE mo_integrals_map mo_bielec_integral_schwartz mo_bielec_integral_jj mo_bielec_integral_jj_anti + FREE mo_integrals_map mo_bielec_integral_jj mo_bielec_integral_jj_anti FREE mo_bielec_integral_jj_exchange mo_bielec_integrals_in_map @@ -1478,7 +1462,7 @@ end subroutine provide_all_mo_integrals implicit none - provide mo_integrals_map mo_bielec_integral_schwartz mo_bielec_integral_jj mo_bielec_integral_jj_anti + provide mo_integrals_map mo_bielec_integral_jj mo_bielec_integral_jj_anti provide mo_bielec_integral_jj_exchange mo_bielec_integrals_in_map end