added check_sym

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:

src/cisd/cisd.irp.f

@@ -53,7 +53,11 @@ subroutine run
   implicit none
   integer :: i
 
+  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:'

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
 

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
+

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)

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