subroutine calc_determ_hamil_sparse_hphf(rep)
type(core_space_t), intent(inout) :: rep
integer :: i, j, row_size, counter, ierr
integer :: nI(nel), nJ(nel)
integer(n_int), allocatable, dimension(:, :) :: temp_store
integer, allocatable :: temp_store_nI(:, :)
integer(TagIntType) :: HRTag, TempStoreTag
HElement_t(dp), allocatable, dimension(:) :: hamiltonian_row
HElement_t(dp) :: HOffDiag
character(len=*), parameter :: t_r = "calc_determ_hamil_sparse_hphf"
integer(n_int) :: tmp(0:NIfD)
integer :: IC
logical :: CS_I
logical, allocatable :: cs(:)
allocate(rep%sparse_core_ham(rep%determ_sizes(iProcIndex)), stat=ierr)
allocate(rep%SparseCoreHamilTags(2, rep%determ_sizes(iProcIndex)))
allocate(hamiltonian_row(rep%determ_space_size), stat=ierr)
call LogMemAlloc('hamiltonian_row', int(rep%determ_space_size), 8, t_r, HRTag, ierr)
allocate(rep%core_ham_diag(rep%determ_sizes(iProcIndex)), stat=ierr)
allocate(temp_store(0:NIfTot, rep%determ_space_size), stat=ierr)
call LogMemAlloc('temp_store', rep%determ_space_size * (NIfTot + 1), 8, t_r, TempStoreTag, ierr)
allocate(temp_store_nI(nel, rep%determ_space_size), stat=ierr)
allocate(cs(rep%determ_space_size), stat=ierr)
! Stick together the deterministic states from all processors, on
! all processors.
! n.b. Explicitly use 0:NIfTot, as NIfTot may not equal NIfBCast
call MPIAllGatherV(SpawnedParts(0:NIfTot, 1:rep%determ_sizes(iProcIndex)), &
temp_store, rep%determ_sizes, rep%determ_displs)
do i = 1, rep%determ_space_size
call decode_bit_det(temp_store_nI(:, i), temp_store(:, i))
cs(i) = TestClosedShellDet(temp_store(:, i))
end do
! Loop over all deterministic states on this processor.
do i = 1, rep%determ_sizes(iProcIndex)
!call decode_bit_det(nI, SpawnedParts(:, i))
nI = temp_store_nI(:, i + rep%determ_displs(iProcIndex))
row_size = 0
hamiltonian_row = 0.0_dp
CS_I = cs(i + rep%determ_displs(iProcIndex))
! Loop over all deterministic states.
do j = 1, rep%determ_space_size
!call decode_bit_det(nJ, temp_store(:,j))
nJ = temp_store_nI(:, j)
! If on the diagonal of the Hamiltonian.
if (j == i + rep%determ_displs(iProcIndex)) then
hamiltonian_row(j) = get_diagonal_matel(nI, SpawnedParts(:, i)) - Hii
rep%core_ham_diag(i) = hamiltonian_row(j)
! We calculate and store the diagonal matrix element at
! this point for later access.
if (.not. tReadPops) then
call set_det_diagH(i, real(hamiltonian_row(j), dp))
HOffDiag = get_off_diagonal_matel(nI, SpawnedParts(:, i))
call set_det_offdiagH(i, HOffDiag)
end if
! Always include the diagonal elements.
row_size = row_size + 1
else
tmp = ieor(SpawnedParts(0:NIfD, i), temp_store(0:NIfD, j))
tmp = iand(SpawnedParts(0:NIfD, i), tmp)
IC = CountBits(tmp, NIfD)
if (IC <= maxExcit .or. ((.not. CS_I) .and. (.not. cs(j)))) then
if(t_evolve_adjoint(rep%first_run())) then
hamiltonian_row(j) = hphf_off_diag_helement_opt(nJ, temp_store(:, j), SpawnedParts(:, i), IC, cs(j), CS_I)
else
hamiltonian_row(j) = hphf_off_diag_helement_opt(nI, SpawnedParts(:, i), temp_store(:, j), IC, CS_I, cs(j))
endif
if (abs(hamiltonian_row(j)) > 0.0_dp) row_size = row_size + 1
end if
end if
end do
! Now we know the number of non-zero elements in this row of the Hamiltonian, so allocate it.
call allocate_sparse_ham_row(rep%sparse_core_ham, i, row_size, "sparse_core_ham", rep%SparseCoreHamilTags(:, i))
rep%sparse_core_ham(i)%elements = 0.0_dp
rep%sparse_core_ham(i)%positions = 0
rep%sparse_core_ham(i)%num_elements = row_size
counter = 1
do j = 1, rep%determ_space_size
! If non-zero or a diagonal element.
if (abs(hamiltonian_row(j)) > 0.0_dp .or. (j == i + rep%determ_displs(iProcIndex))) then
rep%sparse_core_ham(i)%positions(counter) = j
rep%sparse_core_ham(i)%elements(counter) = hamiltonian_row(j)
counter = counter + 1
end if
if (counter == row_size + 1) exit
end do
end do
! Don't time the mpi_barrier call, because if we did then we wouldn't
! be able separate out some of the core Hamiltonian creation time from
! the MPIBarrier calls in the main loop.
call MPIBarrier(ierr, tTimeIn=.false.)
deallocate(temp_store, stat=ierr)
call LogMemDealloc(t_r, TempStoreTag, ierr)
deallocate(hamiltonian_row, stat=ierr)
call LogMemDealloc(t_r, HRTag, ierr)
deallocate(temp_store_nI, stat=ierr)
deallocate(cs, stat=ierr)
end subroutine calc_determ_hamil_sparse_hphf