From 0b3c2804a231af336b10a7e353e8fdaeb1fd3731 Mon Sep 17 00:00:00 2001 From: eginer Date: Fri, 26 Apr 2019 17:31:15 +0200 Subject: [PATCH] added check_sym --- scripts/generate_h_apply.py | 16 +++++++++++-- src/cisd/cisd.irp.f | 6 ++++- src/cisd/h_apply.irp.f | 4 ++++ src/determinants/EZFIO.cfg | 13 +++++++++++ src/determinants/h_apply.template.f | 5 +++++ src/determinants/slater_rules.irp.f | 35 +++++++++++++++++++++++++++++ 6 files changed, 76 insertions(+), 3 deletions(-) diff --git a/scripts/generate_h_apply.py b/scripts/generate_h_apply.py index b3e88fdb..64f05223 100644 --- a/scripts/generate_h_apply.py +++ b/scripts/generate_h_apply.py @@ -5,6 +5,8 @@ import os keywords = """ check_double_excitation copy_buffer +filter_only_connected_to_hf_single +filter_only_connected_to_hf_double declarations decls_main deinit_thread @@ -256,14 +258,14 @@ class H_apply(object): def filter_only_1h1p(self): self["filter_only_1h1p_single"] = """ -! ! DIR$ FORCEINLINE if (is_a_1h1p(hole).eqv..False.) cycle """ self["filter_only_1h1p_double"] = """ -! ! DIR$ FORCEINLINE if (is_a_1h1p(key).eqv..False.) cycle """ + + def filter_only_2h2p(self): self["filter_only_2h2p_single"] = """ ! ! DIR$ FORCEINLINE @@ -294,6 +296,16 @@ class H_apply(object): if (is_a_two_holes_two_particles(hole)) cycle """ + def filter_only_connected_to_hf(self): + self["filter_only_connected_to_hf_single"] = """ + call connected_to_hf(hole,yes_no) + if (yes_no.eqv..False.) cycle + """ + self["filter_only_connected_to_hf_double"] = """ + call connected_to_hf(key,yes_no) + if (yes_no.eqv..False.) cycle + """ + def set_perturbation(self,pert): if self.perturbation is not None: diff --git a/src/cisd/cisd.irp.f b/src/cisd/cisd.irp.f index 9631f971..65f943d3 100644 --- a/src/cisd/cisd.irp.f +++ b/src/cisd/cisd.irp.f @@ -53,7 +53,11 @@ subroutine run implicit none integer :: i - call H_apply_cisd + if(pseudo_sym)then + call H_apply_cisd_sym + else + call H_apply_cisd + endif print *, 'N_det = ', N_det print*,'******************************' print *, 'Energies of the states:' diff --git a/src/cisd/h_apply.irp.f b/src/cisd/h_apply.irp.f index 1c864607..b18f6373 100644 --- a/src/cisd/h_apply.irp.f +++ b/src/cisd/h_apply.irp.f @@ -5,5 +5,9 @@ BEGIN_SHELL [ /usr/bin/env python2 ] from generate_h_apply import H_apply H = H_apply("cisd",do_double_exc=True) print H + +H = H_apply("cisd_sym",do_double_exc=True) +H.filter_only_connected_to_hf() +print H END_SHELL diff --git a/src/determinants/EZFIO.cfg b/src/determinants/EZFIO.cfg index 12a2c2d1..0b289e99 100644 --- a/src/determinants/EZFIO.cfg +++ b/src/determinants/EZFIO.cfg @@ -89,3 +89,16 @@ doc: Weight of the states in state-average calculations. interface: ezfio size: (determinants.n_states) + +[thresh_sym] +type: Threshold +doc: Thresholds to check if a determinant is connected with HF +interface: ezfio,provider,ocaml +default: 1.e-15 + +[pseudo_sym] +type: logical +doc: If |true|, discard any Slater determinants with an interaction smaller than thresh_sym with HF. +interface: ezfio,provider,ocaml +default: False + diff --git a/src/determinants/h_apply.template.f b/src/determinants/h_apply.template.f index 169bd47e..f16d2e7f 100644 --- a/src/determinants/h_apply.template.f +++ b/src/determinants/h_apply.template.f @@ -150,6 +150,7 @@ subroutine $subroutine_diexcOrg(key_in,key_mask,hole_1,particl_1,hole_2, particl logical :: is_a_2h1p logical :: is_a_2h logical :: b_cycle + logical :: yes_no check_double_excitation = .True. iproc = iproc_in @@ -284,6 +285,7 @@ subroutine $subroutine_diexcOrg(key_in,key_mask,hole_1,particl_1,hole_2, particl $only_1h_double $only_1p_double $only_2h1p_double + $filter_only_connected_to_hf_double key_idx += 1 do k=1,N_int keys_out(k,1,key_idx) = key(k,1) @@ -339,6 +341,7 @@ subroutine $subroutine_diexcOrg(key_in,key_mask,hole_1,particl_1,hole_2, particl $only_1h_double $only_1p_double $only_2h1p_double + $filter_only_connected_to_hf_double key_idx += 1 do k=1,N_int keys_out(k,1,key_idx) = key(k,1) @@ -412,6 +415,7 @@ subroutine $subroutine_monoexc(key_in, hole_1,particl_1,fock_diag_tmp,i_generato logical :: is_a_1h logical :: is_a_1p logical :: is_a_2p + logical :: yes_no do k=1,N_int key_mask(k,1) = 0_bit_kind @@ -493,6 +497,7 @@ subroutine $subroutine_monoexc(key_in, hole_1,particl_1,fock_diag_tmp,i_generato $filter_only_1h1p_single $filter_only_1h2p_single $filter_only_2h2p_single + $filter_only_connected_to_hf_single key_idx += 1 do k=1,N_int keys_out(k,1,key_idx) = hole(k,1) diff --git a/src/determinants/slater_rules.irp.f b/src/determinants/slater_rules.irp.f index 3f0f9335..6b164816 100644 --- a/src/determinants/slater_rules.irp.f +++ b/src/determinants/slater_rules.irp.f @@ -2257,3 +2257,38 @@ subroutine i_H_j_double_alpha_beta(key_i,key_j,Nint,hij) end +subroutine connected_to_hf(key_i,yes_no) + implicit none + use bitmasks + integer(bit_kind), intent(in) :: key_i(N_int,2) + logical , intent(out) :: yes_no + double precision :: hij,phase + integer :: exc(0:2,2,2) + integer :: degree + integer :: m,p + yes_no = .True. + call get_excitation_degree(ref_bitmask,key_i,degree,N_int) + if(degree == 2)then + call i_H_j(ref_bitmask,key_i,N_int,hij) + if(dabs(hij) .lt. thresh_sym)then + yes_no = .False. + endif + else if(degree == 1)then + call get_single_excitation(ref_bitmask,key_i,exc,phase,N_int) + ! Single alpha + if (exc(0,1,1) == 1) then + m = exc(1,1,1) + p = exc(1,2,1) + ! Single beta + else + m = exc(1,1,2) + p = exc(1,2,2) + endif + hij = mo_one_e_integrals(m,p) + if(dabs(hij) .lt. thresh_sym)then + yes_no = .False. + endif + else if(degree == 0)then + yes_no = .True. + endif +end