Fixed multi-state

plugins/Full_CI_ZMQ/pt2_stoch_routines.irp.f

@@ -25,8 +25,8 @@ subroutine ZMQ_pt2(E, pt2,relative_error, absolute_error, error)
   
   double precision               :: sumabove(comb_teeth), sum2above(comb_teeth), Nabove(comb_teeth)
   double precision, external     :: omp_get_wtime
+  double precision               :: state_average_weight_save(N_states), w(N_states)
   double precision               :: time
-  double precision               :: w(N_states)
   integer(ZMQ_PTR), external     :: new_zmq_to_qp_run_socket
   
   if (N_det < max(10,N_states)) then
@@ -35,11 +35,12 @@ subroutine ZMQ_pt2(E, pt2,relative_error, absolute_error, error)
     error(:) = 0.d0
   else
     
+    state_average_weight_save(:) = state_average_weight(:)
     do pt2_stoch_istate=1,N_states
       SOFT_TOUCH pt2_stoch_istate
-      w(:) = 0.d0
+      state_average_weight(:) = 0.d0
-      w(pt2_stoch_istate) = 1.d0
+      state_average_weight(pt2_stoch_istate) = 1.d0
-      call update_psi_average_norm_contrib(w)
+      TOUCH state_average_weight 
       
       allocate(pt2_detail(N_states,N_det_generators+1), comb(N_det_generators), computed(N_det_generators), tbc(0:size_tbc))
       sumabove = 0d0
@@ -141,7 +142,9 @@ subroutine ZMQ_pt2(E, pt2,relative_error, absolute_error, error)
       
       deallocate(pt2_detail, comb, computed, tbc)
     enddo
-    FREE psi_average_norm_contrib pt2_stoch_istate
+    FREE pt2_stoch_istate 
+    state_average_weight(:) = state_average_weight_save(:)
+    TOUCH state_average_weight
   endif
   do k=N_det+1,N_states
     pt2(k) = 0.d0

plugins/Generators_full/generators.irp.f

@@ -30,15 +30,8 @@ END_PROVIDER
  ! Hartree-Fock determinant
  END_DOC
  integer                        :: i, k
- psi_coef_generators = 0.d0
+ psi_det_generators(1:N_int,1:2,1:N_det) = psi_det_sorted(1:N_int,1:2,1:N_det)
- psi_det_generators  = 0_bit_kind
+ psi_coef_generators(1:N_det,1:N_states) = psi_coef_sorted(1:N_det,1:N_states)
- do i=1,N_det_generators
-   do k=1,N_int
-     psi_det_generators(k,1,i) = psi_det_sorted(k,1,i)
-     psi_det_generators(k,2,i) = psi_det_sorted(k,2,i)
-   enddo
-   psi_coef_generators(i,:) = psi_coef_sorted(i,:)
- enddo
 
 END_PROVIDER
 

src/Determinants/density_matrix.irp.f

@@ -368,13 +368,13 @@ BEGIN_PROVIDER [ double precision, state_average_weight, (N_states) ]
  END_DOC
  logical :: exists
 
- state_average_weight = 1.d0
+ state_average_weight(:) = 1.d0
  call ezfio_has_determinants_state_average_weight(exists)
  if (exists) then
   call ezfio_get_determinants_state_average_weight(state_average_weight)
  endif
- state_average_weight = state_average_weight+1.d-31
+ state_average_weight(:) = state_average_weight(:)+1.d-31
- state_average_weight = state_average_weight/(sum(state_average_weight))
+ state_average_weight(:) = state_average_weight(:)/(sum(state_average_weight(:)))
 END_PROVIDER
 
 

src/Determinants/determinants.irp.f

@@ -225,34 +225,6 @@ BEGIN_PROVIDER [ double precision, psi_coef, (psi_det_size,N_states) ]
   
 END_PROVIDER
 
-subroutine update_psi_average_norm_contrib(w)
-  implicit none
-  BEGIN_DOC
-  ! Compute psi_average_norm_contrib for different state average weights w(:)
-  END_DOC
-  double precision, intent(in)   :: w(N_states)
-  double precision               :: w0(N_states), f
-  w0(:) = w(:)/sum(w(:))
-  
-  integer                        :: i,j,k
-  do i=1,N_det
-    psi_average_norm_contrib(i) = psi_coef(i,1)*psi_coef(i,1)*w(1)
-  enddo
-  do k=2,N_states
-    do i=1,N_det
-      psi_average_norm_contrib(i) = psi_average_norm_contrib(i) +    &
-          psi_coef(i,k)*psi_coef(i,k)*w(k)
-    enddo
-  enddo
-  f = 1.d0/sum(psi_average_norm_contrib(1:N_det))
-  do i=1,N_det
-    psi_average_norm_contrib(i) = psi_average_norm_contrib(i)*f
-  enddo
-  SOFT_TOUCH psi_average_norm_contrib
-  
-end subroutine
-
-
 BEGIN_PROVIDER [ double precision, psi_average_norm_contrib, (psi_det_size) ]
   implicit none
   BEGIN_DOC
@@ -260,14 +232,12 @@ BEGIN_PROVIDER [ double precision, psi_average_norm_contrib, (psi_det_size) ]
   END_DOC
   integer                        :: i,j,k
   double precision               :: f
-  f = 1.d0/dble(N_states)
+
-  do i=1,N_det
+  psi_average_norm_contrib(:) = 0.d0
-    psi_average_norm_contrib(i) = psi_coef(i,1)*psi_coef(i,1)*f
+  do k=1,N_states
-  enddo
-  do k=2,N_states
     do i=1,N_det
       psi_average_norm_contrib(i) = psi_average_norm_contrib(i) +    &
-          psi_coef(i,k)*psi_coef(i,k)*f
+          psi_coef(i,k)*psi_coef(i,k)*state_average_weight(k)
     enddo
   enddo
   f = 1.d0/sum(psi_average_norm_contrib(1:N_det))