subroutine calc_approx_hamil_sparse_hphf(rep)
type(core_space_t), intent(in) :: rep
integer :: i, j, row_size, counter, ierr
integer :: nI(nel), nJ(nel)
integer, allocatable :: temp_store_nI(:, :)
HElement_t(dp), allocatable, dimension(:) :: hamiltonian_row
HElement_t(dp) :: HOffDiag
character(len=*), parameter :: t_r = "calc_approx_hamil_sparse_hphf"
integer(n_int) :: tmp(0:NIfD)
integer :: IC
logical :: CS_I, var_state_i, var_state_j
logical, allocatable :: cs(:)
allocate(approx_ham(rep%determ_sizes(iProcIndex)), stat=ierr)
allocate(hamiltonian_row(rep%determ_space_size), stat=ierr)
allocate(temp_store_nI(nel, rep%determ_space_size), stat=ierr)
allocate(cs(rep%determ_space_size), stat=ierr)
do i = 1, rep%determ_space_size
call decode_bit_det(temp_store_nI(:, i), rep%core_space(:, i))
cs(i) = TestClosedShellDet(rep%core_space(:, 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))
var_state_i = is_var_state(rep%core_space(:, i + rep%determ_displs(iProcIndex)), nI)
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, core_space(:,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, rep%core_space(:, i + rep%determ_displs(iProcIndex))) - Hii
!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, rep%core_space(:, i + rep%determ_displs(iProcIndex) ) )
call set_det_offdiagH(i, HOffDiag)
end if
! Always include the diagonal elements.
row_size = row_size + 1
else
var_state_j = is_var_state(rep%core_space(:, j), nJ)
! Only add a matrix element if both states are variational
! states, or if one of them is (the rectangular portion of
! H connected var_space to the space of connections to it).
if (var_state_i .or. var_state_j) then
tmp = ieor(rep%core_space(0:NIfD, i + rep%determ_displs(iProcIndex)), rep%core_space(0:NIfD, j))
tmp = iand(rep%core_space(0:NIfD, i + rep%determ_displs(iProcIndex)), 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, &
rep%core_space(:, j), &
rep%core_space(:, i + rep%determ_displs(iProcIndex)), &
IC, cs(j), CS_I)
else
hamiltonian_row(j) = hphf_off_diag_helement_opt(nI, &
rep%core_space(:, i + rep%determ_displs(iProcIndex)), &
rep%core_space(:, 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 if
end do
! Now we know the number of non-zero elements in this row of the Hamiltonian, so allocate it.
allocate(approx_ham(i)%elements(row_size), stat=ierr)
allocate(approx_ham(i)%positions(row_size), stat=ierr)
approx_ham(i)%elements = 0.0_dp
approx_ham(i)%positions = 0
approx_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
approx_ham(i)%positions(counter) = j
approx_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(hamiltonian_row, stat=ierr)
deallocate(temp_store_nI, stat=ierr)
deallocate(cs, stat=ierr)
end subroutine calc_approx_hamil_sparse_hphf