Fixed Full CAS-SD

2024-06-26 15:12:14 +02:00 · 2015-05-06 17:01:45 +02:00 · 2015-05-06 17:01:45 +02:00 · d6ec67e840
commit d6ec67e840
parent bfddb4f9d9
7 changed files with 92 additions and 68 deletions
--- a/src/CAS_SD/cas_sd.irp.f
+++ b/src/CAS_SD/cas_sd.irp.f
@ -69,18 +69,4 @@ program full_ci
    print *,  ''
  enddo
  print *,  'Max excitation degree in the CAS :', exc_max
-  do i=1,N_det
-    in_cas = .False.
-    degree_min = 1000
-    do k=1,N_det_generators
-      call get_excitation_degree(psi_det_generators(1,1,k),psi_det(1,1,i),degree,N_int)
-      degree_min = min(degree_min,degree)
-    enddo
-    if (degree_min > 2) then
-        print *,  'Error : This is not a CAS-SD : '
-        print *,  'Excited determinant:', degree_min
-        call debug_det(psi_det(1,1,k),N_int)
-        stop
-    endif
-  enddo
 end
--- a/src/Determinants/H_apply.irp.f
+++ b/src/Determinants/H_apply.irp.f
@ -53,13 +53,17 @@ subroutine resize_H_apply_buffer(new_size,iproc)
  double precision,  pointer     :: buffer_e2(:,:)
  integer                        :: i,j,k
  integer                        :: Ndet
+
+  BEGIN_DOC
+! Resizes the H_apply buffer of proc iproc. The buffer lock should
+! be set before calling this function.
+  END_DOC
  PROVIDE H_apply_buffer_allocated
  
  ASSERT (new_size > 0)
  ASSERT (iproc >= 0)
  ASSERT (iproc < nproc)
  
-  call omp_set_lock(H_apply_buffer_lock(1,iproc))
  allocate ( buffer_det(N_int,2,new_size),                           &
      buffer_coef(new_size,N_states),                                &
      buffer_e2(new_size,N_states) )
@ -93,7 +97,6 @@ subroutine resize_H_apply_buffer(new_size,iproc)
  
  H_apply_buffer(iproc)%sze = new_size
  H_apply_buffer(iproc)%N_det = min(new_size,H_apply_buffer(iproc)%N_det)
-  call omp_unset_lock(H_apply_buffer_lock(1,iproc))
  
 end

@ -101,8 +104,7 @@ subroutine copy_H_apply_buffer_to_wf
  use omp_lib
  implicit none
  BEGIN_DOC
-! Copies the H_apply buffer to psi_coef. You need to touch psi_det, psi_coef and N_det
-! after calling this function.
+! Copies the H_apply buffer to psi_coef.
 ! After calling this subroutine, N_det, psi_det and psi_coef need to be touched
  END_DOC
  integer(bit_kind), allocatable :: buffer_det(:,:,:)
@ -181,42 +183,77 @@ subroutine copy_H_apply_buffer_to_wf
  call normalize(psi_coef,N_det)
  SOFT_TOUCH N_det psi_det psi_coef
  
-  call debug_unicity_of_determinants
+  logical :: found_duplicates
+  call remove_duplicates_in_psi_det(found_duplicates)
 end

-subroutine debug_unicity_of_determinants
+subroutine remove_duplicates_in_psi_det(found_duplicates)
  implicit none
+  logical, intent(out) :: found_duplicates
  BEGIN_DOC
-! This subroutine checks that there are no repetitions in the wave function
+! Removes duplicate determinants in the wave function.
  END_DOC
-  logical :: same, failed
-  integer :: i,k
-  print *,  "======= DEBUG UNICITY ========="
-  failed = .False.
-  do i=2,N_det
-    same = .True.
-    do k=1,N_int
-      if ( psi_det_sorted_bit(k,1,i) /= psi_det_sorted_bit(k,1,i-1) ) then
-        same = .False.
-        exit
+  integer                        :: i,j,k
+  integer(bit_kind), allocatable :: bit_tmp(:)
+  logical,allocatable            :: duplicate(:)
+
+  allocate (duplicate(N_det), bit_tmp(N_det))
+
+  do i=1,N_det
+    integer, external            :: det_search_key
+    !$DIR FORCEINLINE
+    bit_tmp(i) = det_search_key(psi_det_sorted_bit(1,1,i),N_int)
+    duplicate(i) = .False.
+  enddo
+
+  do i=1,N_det-1
+    if (duplicate(i)) then
+      cycle
+    endif
+    j = i+1
+    do while (bit_tmp(j)==bit_tmp(i))
+      if (duplicate(j)) then
+        j += 1
+        cycle
      endif
-      if ( psi_det_sorted_bit(k,2,i) /= psi_det_sorted_bit(k,2,i-1) ) then
-        same = .False.
+      duplicate(j) = .True.
+      do k=1,N_int
+        if ( (psi_det_sorted_bit(k,1,i) /= psi_det_sorted_bit(k,1,j) ) &
+        .or. (psi_det_sorted_bit(k,2,i) /= psi_det_sorted_bit(k,2,j) ) ) then
+          duplicate(j) = .False.
+          exit
+        endif
+      enddo
+      j += 1
+      if (j > N_det) then
        exit
      endif
    enddo
-    if (same) then
-      failed = .True.
-      call debug_det(psi_det_sorted_bit(1,1,i))
+  enddo
+
+  found_duplicates = .False.
+  do i=1,N_det
+    if (duplicate(i)) then
+      found_duplicates = .True.
+      exit
    endif
  enddo

-  if (failed) then
-    print *,  '======= Determinants not unique : Failed ! ========='
-    stop
-  else
-    print *,  '======= Determinants are unique : OK ! ========='
+  if (found_duplicates) then
+    call write_bool(output_determinants,found_duplicates,'Found duplicate determinants')
+    k=0
+    do i=1,N_det
+      if (.not.duplicate(i)) then
+        k += 1
+        psi_det(:,:,k) = psi_det_sorted_bit (:,:,i)
+        psi_coef(k,:)  = psi_coef_sorted_bit(i,:)
+      endif
+    enddo
+    print *,  N_det,k
+    N_det = k
+    TOUCH N_det psi_det psi_coef
  endif
+  deallocate (duplicate,bit_tmp)
 end

 subroutine fill_H_apply_buffer_no_selection(n_selected,det_buffer,Nint,iproc)
@ -231,11 +268,11 @@ subroutine fill_H_apply_buffer_no_selection(n_selected,det_buffer,Nint,iproc)
  integer                        :: i,j,k
  integer                        :: new_size
  PROVIDE H_apply_buffer_allocated
+  call omp_set_lock(H_apply_buffer_lock(1,iproc))
  new_size = H_apply_buffer(iproc)%N_det + n_selected
  if (new_size > H_apply_buffer(iproc)%sze) then
    call resize_h_apply_buffer(max(2*H_apply_buffer(iproc)%sze,new_size),iproc)
  endif
-  call omp_set_lock(H_apply_buffer_lock(1,iproc))
  do i=1,H_apply_buffer(iproc)%N_det
    ASSERT (sum(popcnt(H_apply_buffer(iproc)%det(:,1,i)) )== elec_alpha_num)
    ASSERT (sum(popcnt(H_apply_buffer(iproc)%det(:,2,i))) == elec_beta_num)
@ -250,7 +287,7 @@ subroutine fill_H_apply_buffer_no_selection(n_selected,det_buffer,Nint,iproc)
  enddo
  do j=1,N_states
    do i=1,N_selected
-      H_apply_buffer(iproc)%coef(i,j) = 0.d0
+      H_apply_buffer(iproc)%coef(i+H_apply_buffer(iproc)%N_det,j) = 0.d0
    enddo
  enddo
  H_apply_buffer(iproc)%N_det = new_size
--- a/src/Determinants/README.rst
+++ b/src/Determinants/README.rst
@ -187,10 +187,10 @@ Documentation
 `det_svd <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/det_svd.irp.f#L1>`_
  Computes the SVD of the Alpha x Beta determinant coefficient matrix

-`filter_3_highest_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L426>`_
+`filter_3_highest_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L423>`_
  Returns a determinant with only the 3 highest electrons

-`int_of_3_highest_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L391>`_
+`int_of_3_highest_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L388>`_
  Returns an integer*8 as :
  .br
  |_<--- 21 bits ---><--- 21 bits ---><--- 21 bits --->|
@ -207,26 +207,26 @@ Documentation
 `n_det <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L3>`_
  Number of determinants in the wave function

-`psi_average_norm_contrib <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L276>`_
+`psi_average_norm_contrib <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L273>`_
  Contribution of determinants to the state-averaged density

-`psi_average_norm_contrib_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L306>`_
+`psi_average_norm_contrib_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L303>`_
  Wave function sorted by determinants contribution to the norm (state-averaged)

-`psi_coef <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L230>`_
+`psi_coef <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L227>`_
  The wave function coefficients. Initialized with Hartree-Fock if the EZFIO file
  is empty

-`psi_coef_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L305>`_
+`psi_coef_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L302>`_
  Wave function sorted by determinants contribution to the norm (state-averaged)

-`psi_coef_sorted_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L453>`_
+`psi_coef_sorted_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L450>`_
  Determinants on which we apply <i|H|j>.
  They are sorted by the 3 highest electrons in the alpha part,
  then by the 3 highest electrons in the beta part to accelerate
  the research of connected determinants.

-`psi_coef_sorted_bit <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L336>`_
+`psi_coef_sorted_bit <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L333>`_
  Determinants on which we apply <i|H|psi> for perturbation.
  They are sorted by determinants interpreted as integers. Useful
  to accelerate the search of a random determinant in the wave
@ -239,46 +239,46 @@ Documentation
 `psi_det_size <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L47>`_
  Size of the psi_det/psi_coef arrays

-`psi_det_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L304>`_
+`psi_det_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L301>`_
  Wave function sorted by determinants contribution to the norm (state-averaged)

-`psi_det_sorted_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L452>`_
+`psi_det_sorted_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L449>`_
  Determinants on which we apply <i|H|j>.
  They are sorted by the 3 highest electrons in the alpha part,
  then by the 3 highest electrons in the beta part to accelerate
  the research of connected determinants.

-`psi_det_sorted_bit <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L335>`_
+`psi_det_sorted_bit <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L332>`_
  Determinants on which we apply <i|H|psi> for perturbation.
  They are sorted by determinants interpreted as integers. Useful
  to accelerate the search of a random determinant in the wave
  function.

-`psi_det_sorted_next_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L454>`_
+`psi_det_sorted_next_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L451>`_
  Determinants on which we apply <i|H|j>.
  They are sorted by the 3 highest electrons in the alpha part,
  then by the 3 highest electrons in the beta part to accelerate
  the research of connected determinants.

-`read_dets <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L583>`_
+`read_dets <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L580>`_
  Reads the determinants from the EZFIO file

-`save_wavefunction <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L630>`_
+`save_wavefunction <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L627>`_
  Save the wave function into the EZFIO file

-`save_wavefunction_general <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L649>`_
+`save_wavefunction_general <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L646>`_
  Save the wave function into the EZFIO file

-`save_wavefunction_unsorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L640>`_
+`save_wavefunction_unsorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L637>`_
  Save the wave function into the EZFIO file

-`sort_dets_by_3_highest_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L474>`_
+`sort_dets_by_3_highest_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L471>`_
  Determinants on which we apply <i|H|j>.
  They are sorted by the 3 highest electrons in the alpha part,
  then by the 3 highest electrons in the beta part to accelerate
  the research of connected determinants.

-`sort_dets_by_det_search_key <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L349>`_
+`sort_dets_by_det_search_key <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L346>`_
  Determinants are sorted are sorted according to their det_search_key.
  Useful to accelerate the search of a random determinant in the wave
  function.
@ -316,7 +316,7 @@ Documentation
 `diag_algorithm <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI.irp.f#L1>`_
  Diagonalization algorithm (Davidson or Lapack)

-`diagonalize_ci <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI.irp.f#L96>`_
+`diagonalize_ci <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI.irp.f#L100>`_
  Replace the coefficients of the CI states by the coefficients of the
  eigenstates of the CI matrix

@ -345,7 +345,7 @@ Documentation
 `ci_electronic_energy_mono <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI_mono.irp.f#L1>`_
  Eigenvectors/values of the CI matrix

-`diagonalize_ci_mono <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI_mono.irp.f#L59>`_
+`diagonalize_ci_mono <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI_mono.irp.f#L73>`_
  Replace the coefficients of the CI states by the coefficients of the
  eigenstates of the CI matrix

--- a/src/Determinants/s2.irp.f
+++ b/src/Determinants/s2.irp.f
@ -88,7 +88,7 @@ subroutine get_s2_u0(psi_keys_tmp,psi_coefs_tmp,n,nmax,s2)
 double precision, intent(out) :: s2
 integer :: i,j,l
 double precision :: s2_tmp
- s2 = S_z2_Sz
+ s2 = 0.d0
 !$OMP PARALLEL DO DEFAULT(NONE) &
 !$OMP PRIVATE(i,j,s2_tmp) SHARED(n,psi_coefs_tmp,psi_keys_tmp,N_int) &
 !$OMP REDUCTION(+:s2) SCHEDULE(dynamic)
@ -101,5 +101,6 @@ subroutine get_s2_u0(psi_keys_tmp,psi_coefs_tmp,n,nmax,s2)
   enddo
 enddo
 !$OMP END PARALLEL DO
+ s2 +=  S_z2_Sz
 end

--- a/src/MRCC/README.rst
+++ b/src/MRCC/README.rst
@ -51,7 +51,7 @@ Documentation
 `ci_electronic_energy_dressed <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC/mrcc_utils.irp.f#L78>`_
  Eigenvectors/values of the CI matrix

-`ci_energy_dressed <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC/mrcc_utils.irp.f#L132>`_
+`ci_energy_dressed <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC/mrcc_utils.irp.f#L144>`_
  N_states lowest eigenvalues of the dressed CI matrix

 `delta_ij <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC/mrcc_utils.irp.f#L43>`_
@ -60,7 +60,7 @@ Documentation
 `delta_ij_non_cas <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC/mrcc_utils.irp.f#L34>`_
  Dressing matrix in SD basis

-`diagonalize_ci_dressed <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC/mrcc_utils.irp.f#L147>`_
+`diagonalize_ci_dressed <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC/mrcc_utils.irp.f#L159>`_
  Replace the coefficients of the CI states by the coefficients of the
  eigenstates of the CI matrix

--- a/src/MonoInts/pot_ao_ints.irp.f
+++ b/src/MonoInts/pot_ao_ints.irp.f
@ -10,7 +10,7 @@
 integer           :: i,j,k,l,n_pt_in,m
 double precision  ::overlap_x,overlap_y,overlap_z,overlap,dx,NAI_pol_mult

-  if (do_pseudo.eqv..TRUE.) then
+  if (do_pseudo) then
    ao_nucl_elec_integral = ao_nucl_elec_integral_pseudo
  else
    ao_nucl_elec_integral = 0.d0
--- a/src/Perturbation/selection.irp.f
+++ b/src/Perturbation/selection.irp.f
@ -19,6 +19,7 @@ subroutine fill_H_apply_buffer_selection(n_selected,det_buffer,e_2_pert_buffer,c
  ASSERT (Nint > 0)
  ASSERT (N_int == N_int)
  ASSERT (N_selected >= 0)
+  call omp_set_lock(H_apply_buffer_lock(1,iproc))
  smax = selection_criterion
  smin = selection_criterion_min
  new_size = H_apply_buffer(iproc)%N_det + n_selected
@ -26,7 +27,6 @@ subroutine fill_H_apply_buffer_selection(n_selected,det_buffer,e_2_pert_buffer,c
  if (new_size > h_apply_buffer(iproc)%sze) then
    call resize_h_apply_buffer(max(h_apply_buffer(iproc)%sze*2,new_size),iproc)
  endif
-  call omp_set_lock(H_apply_buffer_lock(1,iproc))
  do i=1,H_apply_buffer(iproc)%N_det
    ASSERT (sum(popcnt(h_apply_buffer(iproc)%det(:,1,i)) )== elec_alpha_num)
    ASSERT (sum(popcnt(h_apply_buffer(iproc)%det(:,2,i))) == elec_beta_num)