minilist indices in psi_det

plugins/dress_zmq/alpha_factory.irp.f

@@ -20,6 +20,23 @@ subroutine alpha_callback(delta_ij_loc, i_generator, subset,iproc)
 end subroutine
 
 
+ BEGIN_PROVIDER [ integer, psi_from_sorted, (N_det) ]
+&BEGIN_PROVIDER [ integer, idx_non_ref_from_sorted, (N_det) ]
+  implicit none
+  integer :: i,inpsisor
+  
+  idx_non_ref_from_sorted = 0
+  psi_from_sorted = 0
+
+  do i=1,N_det
+    psi_from_sorted(psi_det_sorted_order(i)) = i
+    inpsisor = psi_det_sorted_order(i)
+    if(inpsisor <= 0) stop "idx_non_ref_from_sorted"
+    idx_non_ref_from_sorted(inpsisor) = idx_non_ref_rev(i)
+  end do
+END_PROVIDER
+
+
 subroutine generate_singles_and_doubles(delta_ij_loc, i_generator, bitmask_index, subset, iproc)
   use bitmasks
   implicit none
@@ -371,17 +388,17 @@ subroutine generate_singles_and_doubles(delta_ij_loc, i_generator, bitmask_index
 end subroutine
 
 
-subroutine alpha_callback_mask(delta_ij_loc, sp, mask, bannedOrb, banned, indexes, indexes_end, abuf, siz, iproc)
+subroutine alpha_callback_mask(delta_ij_loc, sp, mask, bannedOrb, banned, indexes, indexes_end, rabuf, siz, iproc)
   use bitmasks
   implicit none
 
   double precision,intent(inout) :: delta_ij_loc(N_states,N_det,2) 
   integer, intent(in) :: sp, indexes(0:mo_tot_num, 0:mo_tot_num), siz, iproc
-  integer, intent(in) :: indexes_end(0:mo_tot_num, 0:mo_tot_num), abuf(*)
+  integer, intent(in) :: indexes_end(0:mo_tot_num, 0:mo_tot_num), rabuf(*)
   logical, intent(in) :: bannedOrb(mo_tot_num,2), banned(mo_tot_num, mo_tot_num)
   integer(bit_kind), intent(in) :: mask(N_int, 2)
   integer(bit_kind) :: alpha(N_int, 2)
-  integer, allocatable :: labuf(:)
+  integer, allocatable :: labuf(:), abuf(:)
   logical :: ok
   integer :: i,j,k,s,st1,st2,st3,st4
   integer :: lindex(mo_tot_num,2), lindex_end(mo_tot_num, 2)
@@ -390,14 +407,19 @@ subroutine alpha_callback_mask(delta_ij_loc, sp, mask, bannedOrb, banned, indexe
   integer(bit_kind), allocatable :: det_minilist(:,:,:)
 
 
-  allocate(labuf(N_det), putten(N_det), det_minilist(N_int, 2, N_det))
+  allocate(abuf(siz), labuf(N_det), putten(N_det), det_minilist(N_int, 2, N_det))
+  
+  do i=1,siz
+    abuf(i) = psi_from_sorted(rabuf(i))
+  end do
+
   putten = .false.
 
   st1 = indexes_end(0,0)-1 !!
   if(st1 > 0) then
     labuf(:st1) = abuf(:st1)
     do i=1,st1
-      det_minilist(:,:,i) = psi_det_sorted(:,:,labuf(i))
+      det_minilist(:,:,i) = psi_det(:,:,labuf(i))
     end do
   end if
   st1 += 1
@@ -421,66 +443,66 @@ subroutine alpha_callback_mask(delta_ij_loc, sp, mask, bannedOrb, banned, indexe
     lindex_end(:,1) = indexes_end(1:, 0)-1
   end if
 
-    do i=1,mo_tot_num
+  do i=1,mo_tot_num
-      if(bannedOrb(i,s1)) cycle
+    if(bannedOrb(i,s1)) cycle
-      if(lindex(i,s1) /= 0) then
+    if(lindex(i,s1) /= 0) then
-        st2 = st1 + 1 + lindex_end(i,s1)-lindex(i,s1)
+      st2 = st1 + 1 + lindex_end(i,s1)-lindex(i,s1)
-        labuf(st1:st2-1) = abuf(lindex(i,s1):lindex_end(i,s1))
+      labuf(st1:st2-1) = abuf(lindex(i,s1):lindex_end(i,s1))
-        do j=st1,st2-1
+      do j=st1,st2-1
-          putten(labuf(j)) = .true.
+        putten(labuf(j)) = .true.
-          det_minilist(:,:,j) = psi_det_sorted(:,:,labuf(j))
+        det_minilist(:,:,j) = psi_det(:,:,labuf(j))
-        end do
-      else
-        st2 = st1
-      end if
-      
-      if(sp == 3) then
-        stamo = 1
-      else
-        stamo = i+1
-      end if
-
-      do j=stamo,mo_tot_num
-        if(bannedOrb(j,s2) .or. banned(i,j)) cycle
-        if(lindex(j,s2) /= 0) then
-          st3 = st2
-          do k=lindex(j,s2), lindex_end(j,s2)
-            if(.not. putten(abuf(k))) then
-              labuf(st3) = abuf(k)
-              det_minilist(:,:,st3) = psi_det_sorted(:,:,abuf(k))
-              st3 += 1
-            end if
-          end do
-        else
-          st3 = st2
-        end if
-
-        if(indexes(i,j) /= 0) then
-          st4 = st3 + 1 + indexes_end(i,j)-indexes(i,j) -1!!
-          labuf(st3:st4-1) = abuf(indexes(i,j):indexes_end(i,j)-1) !!
-          do k=st3, st4-1
-            det_minilist(:,:,k) = psi_det_sorted(:,:,labuf(k))
-          end do
-        else
-          st4 = st3
-        end if
-        !APPLY PART
-        if(st4 > 1) then
-          call apply_particles(mask, s1, i, s2, j, alpha, ok, N_int)
-          !if(.not. ok) stop "non existing alpha......"
-          !print *, "willcall", st4-1, size(labuf)
-          call dress_with_alpha_buffer(delta_ij_loc, labuf, det_minilist, st4-1, alpha, iproc)
-          !call dress_with_alpha_buffer(delta_ij_loc, abuf, siz, alpha, 1)
-        end if
       end do
-      
+    else
-      if(lindex(i,s1) /= 0) then
+      st2 = st1
-        do j=st1,st2-1
+    end if
-          putten(labuf(j)) = .false.
+    
+    if(sp == 3) then
+      stamo = 1
+    else
+      stamo = i+1
+    end if
+
+    do j=stamo,mo_tot_num
+      if(bannedOrb(j,s2) .or. banned(i,j)) cycle
+      if(lindex(j,s2) /= 0) then
+        st3 = st2
+        do k=lindex(j,s2), lindex_end(j,s2)
+          if(.not. putten(abuf(k))) then
+            labuf(st3) = abuf(k)
+            det_minilist(:,:,st3) = psi_det(:,:,abuf(k))
+            st3 += 1
+          end if
         end do
+      else
+        st3 = st2
       end if
 
+      if(indexes(i,j) /= 0) then
+        st4 = st3 + 1 + indexes_end(i,j)-indexes(i,j) -1!!
+        labuf(st3:st4-1) = abuf(indexes(i,j):indexes_end(i,j)-1) !!
+        do k=st3, st4-1
+          det_minilist(:,:,k) = psi_det(:,:,labuf(k))
+        end do
+      else
+        st4 = st3
+      end if
+      !APPLY PART
+      if(st4 > 1) then
+        call apply_particles(mask, s1, i, s2, j, alpha, ok, N_int)
+        !if(.not. ok) stop "non existing alpha......"
+        !print *, "willcall", st4-1, size(labuf)
+        call dress_with_alpha_buffer(delta_ij_loc, labuf, det_minilist, st4-1, alpha, iproc)
+        !call dress_with_alpha_buffer(delta_ij_loc, abuf, siz, alpha, 1)
+      end if
     end do
+    
+    if(lindex(i,s1) /= 0) then
+      do j=st1,st2-1
+        putten(labuf(j)) = .false.
+      end do
+    end if
+
+  end do
 end subroutine
 
 
@@ -543,7 +565,6 @@ subroutine count_pq(mask, sp, det, i_gen, N_sel, bannedOrb, banned, countedGlob,
 
     mobMask(1,1) = iand(negMask(1,1), det(1,1,i))
     mobMask(1,2) = iand(negMask(1,2), det(1,2,i))
-    if(interesting(i) < i_gen) cycle
     nt = popcnt(mobMask(1, 1)) + popcnt(mobMask(1, 2))
 
 
@@ -656,11 +677,11 @@ subroutine splash_pq(mask, sp, det, i_gen, N_sel, bannedOrb, banned, indexes, ab
     !if (interesting(i) < 0) then
     !  stop 'prefetch interesting(i)'
     !endif
+    if(interesting(i) < i_gen) cycle
 
     
     mobMask(1,1) = iand(negMask(1,1), det(1,1,i))
     mobMask(1,2) = iand(negMask(1,2), det(1,2,i))
-    if(interesting(i) < i_gen) cycle
     nt = popcnt(mobMask(1, 1)) + popcnt(mobMask(1, 2))
 
     if(nt > 4) cycle

plugins/mrcc_sto/mrcc_sto.irp.f

@@ -7,28 +7,11 @@ program mrcc_sto
   call dress_zmq()
 end
 
- BEGIN_PROVIDER [ integer, psi_from_sorted, (N_det) ]
-&BEGIN_PROVIDER [ integer, idx_non_ref_from_sorted, (N_det) ]
-  implicit none
-  integer :: i,inpsisor
-  
-  idx_non_ref_from_sorted = 0
-  psi_from_sorted = 0
-
-  do i=1,N_det
-    psi_from_sorted(psi_det_sorted_order(i)) = i
-    inpsisor = psi_det_sorted_order(i)
-    if(inpsisor <= 0) stop "idx_non_ref_from_sorted"
-    idx_non_ref_from_sorted(inpsisor) = idx_non_ref_rev(i)
-  end do
-END_PROVIDER
-
 
  BEGIN_PROVIDER [ double precision, hij_cache_, (N_det,Nproc) ]
 &BEGIN_PROVIDER [ double precision, sij_cache_, (N_det,Nproc) ]
 &BEGIN_PROVIDER [ double precision, dIa_hla_, (N_states,N_det,Nproc) ]
 &BEGIN_PROVIDER [ double precision, dIa_sla_, (N_states,N_det,Nproc) ]
-&BEGIN_PROVIDER [ integer, idx_alpha_, (0:N_det,Nproc) ]
  BEGIN_DOC
  ! temporay arrays for dress_with_alpha_buffer. Avoids realocation.
 END_DOC
@@ -71,7 +54,7 @@ subroutine dress_with_alpha_buffer(delta_ij_loc, minilist, det_minilist, n_minil
   if(n_minilist == 1) return
   
   do i=1,n_minilist
-    if(idx_non_ref_from_sorted(minilist(i)) == 0) return
+    if(idx_non_ref_rev(minilist(i)) == 0) return
   end do
 
   shdress = 0d0
@@ -81,40 +64,20 @@ subroutine dress_with_alpha_buffer(delta_ij_loc, minilist, det_minilist, n_minil
     PROVIDE one_anhil fock_virt_total fock_core_inactive_total one_creat
   endif
   
-  ll_sd = 0
   do l_sd=1,n_minilist
-    ok = .true.
+   !call get_excitation(det_minilist(1,1,l_sd),alpha,exc,degree1,phase,N_int)
-    k_sd = minilist(l_sd)
+   !if(degree1 == 0 .or. degree1 > 2) stop "minilist error"
-    !do i_I=1,N_det_ref
+   !call decode_exc(exc,degree1,h1,p1,h2,p2,s1,s2)
-    !  call get_excitation_degree(psi_det_sorted(1,1,k_sd),psi_ref(1,1,i_I),degree1,N_int)
+   !
-    !  if(degree1 == 0) then
+   !ok = (mo_class(h1)(1:1) == 'A' .or. mo_class(h1)(1:1) == 'I') .and. &
-    !    ok = .false.
+   !     (mo_class(p1)(1:1) == 'A' .or. mo_class(p1)(1:1) == 'V') 
-    !    exit
+   !if(ok .and. degree1 == 2) then
-    !  end if
+   !        ok = (mo_class(h2)(1:1) == 'A' .or. mo_class(h2)(1:1) == 'I') .and. &
-    !end do
+   !             (mo_class(p2)(1:1) == 'A' .or. mo_class(p2)(1:1) == 'V') 
-
-   !if(ok) then
-   !  call get_excitation(psi_det_sorted(1,1,k_sd),alpha,exc,degree1,phase,N_int)
-   !  if(degree1 == 0 .or. degree1 > 2) stop "minilist error"
-   !  call decode_exc(exc,degree1,h1,p1,h2,p2,s1,s2)
-   !  
-   !  ok = (mo_class(h1)(1:1) == 'A' .or. mo_class(h1)(1:1) == 'I') .and. &
-   !       (mo_class(p1)(1:1) == 'A' .or. mo_class(p1)(1:1) == 'V') 
-   !  if(ok .and. degree1 == 2) then
-   !          ok = (mo_class(h2)(1:1) == 'A' .or. mo_class(h2)(1:1) == 'I') .and. &
-   !               (mo_class(p2)(1:1) == 'A' .or. mo_class(p2)(1:1) == 'V') 
-   !  end if
    !end if
-    
+    call i_h_j(alpha,det_minilist(1,1,l_sd),N_int,hij_cache_(l_sd,iproc))
-    if(ok) then
+    call get_s2(alpha,det_minilist(1,1,l_sd),N_int,sij_cache_(l_sd,iproc))
-      ll_sd += 1
-      idx_alpha_(ll_sd,iproc) = k_sd
-      call i_h_j(alpha,psi_det_sorted(1,1,k_sd),N_int,hij_cache_(k_sd,iproc))
-      call get_s2(alpha,psi_det_sorted(1,1,k_sd),N_int,sij_cache_(k_sd,iproc))
-    end if
   enddo
-  if(ll_sd <= 1) return 
-  idx_alpha_(0,iproc) = ll_sd
   
 
   do i_I=1,N_det_ref
@@ -127,62 +90,68 @@ subroutine dress_with_alpha_buffer(delta_ij_loc, minilist, det_minilist, n_minil
       dIa(i_state) = 0.d0
     enddo
 
-    do k_sd=1,idx_alpha_(0,iproc)
+    do k_sd=1,n_minilist
-      call get_excitation_degree(psi_ref(1,1,i_I),psi_det_sorted(1,1,idx_alpha_(k_sd,iproc)),degree,N_int)
+      call get_excitation_degree(psi_ref(1,1,i_I),det_minilist(1,1,k_sd),degree,N_int)
-!      print *, "diden"
       if (degree > 2) then
         cycle
       endif
       
-      call get_excitation(psi_det_sorted(1,1,idx_alpha_(k_sd,iproc)),alpha,exc,degree2,phase,N_int)
+      call get_excitation(det_minilist(1,1,k_sd),alpha,exc,degree2,phase,N_int)
-      !print *, "DEG", degree2
-      call decode_exc(exc,degree2,h1,p1,h2,p2,s1,s2)
-      do k=1,N_int
-        tmp_det(k,1) = psi_ref(k,1,i_I)
-        tmp_det(k,2) = psi_ref(k,2,i_I)
-      enddo
       call apply_excitation(psi_ref(1,1,i_I), exc, tmp_det, ok, N_int)
       
-      ok2 = .false.
+      if((.not. ok) .and. (.not. perturbative_triples)) cycle
-      do i_state=1,N_states
+
-        dIK(i_state) = dij(i_I, idx_non_ref_from_sorted(idx_alpha_(k_sd,iproc)), i_state)
-        if(dIK(i_state) /= 0d0) then
-          ok2 = .true.
-        endif
-      enddo
-      if(.not. ok2) cycle
-      
-      ! <I| \l/ |alpha>
       do i_state=1,N_states
         dka(i_state) = 0.d0
       enddo
+      
+      ok2 = .false.
+      !do i_state=1,N_states
+      !  !if(dka(i_state) == 0) cycle
+      !  dIk(i_state) = dij(i_I, idx_non_ref_rev(minilist(k_sd)), i_state)
+      !  if(dIk(i_state) /= 0d0) then
+      !    ok2 = .true.
+      !  endif
+      !enddo
+      !if(.not. ok2) cycle
 
       if (ok) then
-        do l_sd=k_sd+1,idx_alpha_(0,iproc)
+        phase2 = 0d0
-          call get_excitation_degree(tmp_det,psi_det_sorted(1,1,idx_alpha_(l_sd,iproc)),degree,N_int)
+        do l_sd=k_sd+1,n_minilist
+          call get_excitation_degree(tmp_det,det_minilist(1,1,l_sd),degree,N_int)
           if (degree == 0) then
-            call get_excitation(psi_ref(1,1,i_I),psi_det_sorted(1,1,idx_alpha_(l_sd,iproc)),exc,degree,phase2,N_int)
             do i_state=1,N_states
-              dka(i_state) = dij(i_I, idx_non_ref_from_sorted(idx_alpha_(l_sd,iproc)), i_state) * phase * phase2
+              dIk(i_state) = dij(i_I, idx_non_ref_rev(minilist(k_sd)), i_state)
-            enddo
+              if(dIk(i_state) /= 0d0) then
+                if(phase2 == 0d0) call get_excitation(psi_ref(1,1,i_I),det_minilist(1,1,l_sd),exc,degree,phase2,N_int)
+                dka(i_state) = dij(i_I, idx_non_ref_rev(minilist(l_sd)), i_state) * phase * phase2
+              end if
+            end do
+
+            !call get_excitation(psi_ref(1,1,i_I),det_minilist(1,1,l_sd),exc,degree,phase2,N_int)
+            !do i_state=1,N_states
+            !  if(dIk(i_state) /= 0d0) dka(i_state) = dij(i_I, idx_non_ref_rev(minilist(l_sd)), i_state) * phase * phase2
+            !enddo
             exit
+
           endif
         enddo
       else if (perturbative_triples) then
-          hka = hij_cache_(idx_alpha_(k_sd,iproc),iproc)
+        hka = hij_cache_(k_sd,iproc)
-          if (dabs(hka) > 1.d-12) then
+        if (dabs(hka) > 1.d-12) then
-            call get_delta_e_dyall_general_mp(psi_ref(1,1,i_I),alpha,Delta_E_inv)
+          call get_delta_e_dyall_general_mp(psi_ref(1,1,i_I),alpha,Delta_E_inv)
 
-            do i_state=1,N_states
+          do i_state=1,N_states
-              ASSERT (Delta_E_inv(i_state) < 0.d0)
+            ASSERT (Delta_E_inv(i_state) < 0.d0)
-              dka(i_state) = hka / Delta_E_inv(i_state)
+            dka(i_state) = hka / Delta_E_inv(i_state)
-            enddo
+          enddo
-          endif
+        endif
       endif
-
+      
+           
       if (perturbative_triples.and. (degree2 == 1) ) then
           call i_h_j(psi_ref(1,1,i_I),tmp_det,N_int,hka)
-          hka = hij_cache_(idx_alpha_(k_sd,iproc),iproc) - hka
+          hka = hij_cache_(k_sd,iproc) - hka
           if (dabs(hka) > 1.d-12) then
             call get_delta_e_dyall_general_mp(psi_ref(1,1,i_I),alpha,Delta_E_inv)
             do i_state=1,N_states
@@ -202,32 +171,20 @@ subroutine dress_with_alpha_buffer(delta_ij_loc, minilist, det_minilist, n_minil
     enddo
     if(.not. ok2) cycle
 
-    do l_sd=1,idx_alpha_(0,iproc)
+    do l_sd=1,n_minilist
-      k_sd = idx_alpha_(l_sd,iproc)
+      k_sd = minilist(l_sd)
-      hla = hij_cache_(k_sd,iproc)
+      hla = hij_cache_(l_sd,iproc)
-      sla = sij_cache_(k_sd,iproc)
+      sla = sij_cache_(l_sd,iproc)
       do i_state=1,N_states
-  !     dIa_hla_(i_state,k_sd,iproc) = dIa(i_state) * hla
-  !     dIa_sla_(i_state,k_sd,iproc) = dIa(i_state) * sla
-  !   enddo
-  ! enddo
-  ! do l_sd=1,idx_alpha_(0,iproc)
-  !   do i_state=1,N_states
-        k_sd = idx_alpha_(l_sd,iproc)
-        m_sd = psi_from_sorted(k_sd)
         hdress =  dIa(i_state) * hla * psi_ref_coef(i_I,i_state)
         sdress =  dIa(i_state) * sla * psi_ref_coef(i_I,i_state)
-      ! hdress = dIa_hla_(i_state,k_sd,iproc) * psi_ref_coef(i_I,i_state)
-      ! sdress = dIa_sla_(i_state,k_sd,iproc) * psi_ref_coef(i_I,i_state)
         !$OMP ATOMIC
-        delta_ij_loc(i_state,m_sd,1) += hdress
+        delta_ij_loc(i_state,k_sd,1) += hdress
         !$OMP ATOMIC
-        delta_ij_loc(i_state,m_sd,2) += sdress
+        delta_ij_loc(i_state,k_sd,2) += sdress
-        !print *, "ENDRES" 
       enddo
     enddo
   enddo
-
 end subroutine