Merge branch 'master' into macos

external/qp2-dependencies

@@ -1 +1 @@
-Subproject commit b9e877d9a9444a2aac23ec94b4a174caa4a10dd2
+Subproject commit e0d0e02e9f5ece138d1520106954a881ab0b8db2

src/cipsi/cipsi.irp.f

@@ -130,13 +130,13 @@ subroutine run_cipsi
     if (qp_stop()) exit
   enddo
 
-  if (.not.qp_stop()) then
-    if (N_det < N_det_max) then
-        call diagonalize_CI
-        call save_wavefunction
-        call save_energy(psi_energy_with_nucl_rep, zeros)
-    endif
-
+  ! If stopped because N_det > N_det_max, do an extra iteration to compute the PT2
+  if ((.not.qp_stop()).and.                                          &
+        (N_det > N_det_max) .and.                                    &
+        (maxval(abs(pt2_data % pt2(1:N_states))) > pt2_max) .and.    &
+        (maxval(abs(pt2_data % variance(1:N_states))) > variance_max) .and.&
+        (correlation_energy_ratio <= correlation_energy_ratio_max)   &
+        ) then
     if (do_pt2) then
       call pt2_dealloc(pt2_data)
       call pt2_dealloc(pt2_data_err)

src/cipsi/stochastic_cipsi.irp.f

@@ -119,13 +119,13 @@ subroutine run_stochastic_cipsi
     if (qp_stop()) exit
   enddo
 
-  if (.not.qp_stop()) then
-    if (N_det < N_det_max) then
-        call diagonalize_CI
-        call save_wavefunction
-        call save_energy(psi_energy_with_nucl_rep, zeros)
-    endif
-
+  ! If stopped because N_det > N_det_max, do an extra iteration to compute the PT2
+  if ((.not.qp_stop()).and.                                          &
+        (N_det > N_det_max) .and.                                    &
+        (maxval(abs(pt2_data % pt2(1:N_states))) > pt2_max) .and.    &
+        (maxval(abs(pt2_data % variance(1:N_states))) > variance_max) .and.&
+        (correlation_energy_ratio <= correlation_energy_ratio_max)   &
+        ) then
     call pt2_dealloc(pt2_data)
     call pt2_dealloc(pt2_data_err)
     call pt2_alloc(pt2_data, N_states)

src/davidson/diagonalization_hs2_dressed.irp.f

@@ -465,8 +465,7 @@ subroutine davidson_diag_hjj_sjj(dets_in,u_in,H_jj,s2_out,energies,dim_in,sze,N_
        integer :: lwork, info
        double precision, allocatable :: work(:)
 
-!       y = h
-       y = h_p
+       y = h
        lwork = -1
        allocate(work(1))
        call dsygv(1,'V','U',shift2,y,size(y,1), &

src/two_body_rdm/two_e_dm_mo.irp.f

@@ -16,13 +16,13 @@ BEGIN_PROVIDER [double precision, two_e_dm_mo, (mo_num,mo_num,mo_num,mo_num)]
    two_e_dm_mo = 0.d0
    integer                        :: i,j,k,l,iorb,jorb,korb,lorb,istate
 
-   do l=1,mo_num
+   do l=1,n_core_inact_act_orb
     lorb = list_core_inact_act(l)
-    do k=1,mo_num
+    do k=1,n_core_inact_act_orb
      korb = list_core_inact_act(k)
-     do j=1,mo_num
+     do j=1,n_core_inact_act_orb
       jorb = list_core_inact_act(j)
-      do i=1,mo_num
+      do i=1,n_core_inact_act_orb
         iorb = list_core_inact_act(i)
         two_e_dm_mo(iorb,jorb,korb,lorb) = state_av_full_occ_2_rdm_spin_trace_mo(i,j,k,l)
       enddo