Drain results from local workers and update tracking state.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| type(mbe_context_t), | intent(in) | :: | ctx | |||
| integer(kind=int64), | intent(inout) | :: | worker_fragment_map(:) | |||
| type(calculation_result_t), | intent(inout) | :: | results(:) | |||
| integer(kind=int64), | intent(inout) | :: | results_received | |||
| type(timer_type), | intent(in) | :: | coord_timer |
| Type | Visibility | Attributes | Name | Initial | |||
|---|---|---|---|---|---|---|---|
| logical, | private | :: | has_pending | ||||
| type(MPI_Status), | private | :: | local_status | ||||
| type(request_t), | private | :: | req | ||||
| integer, | private | :: | worker_source |
subroutine handle_local_worker_results(ctx, worker_fragment_map, results, results_received, coord_timer) !! Drain results from local workers and update tracking state. use mqc_many_body_expansion, only: mbe_context_t type(mbe_context_t), intent(in) :: ctx integer(int64), intent(inout) :: worker_fragment_map(:) type(calculation_result_t), intent(inout) :: results(:) integer(int64), intent(inout) :: results_received type(timer_type), intent(in) :: coord_timer type(MPI_Status) :: local_status logical :: has_pending integer :: worker_source type(request_t) :: req do call iprobe(ctx%resources%mpi_comms%node_comm, MPI_ANY_SOURCE, TAG_WORKER_SCALAR_RESULT, has_pending, local_status) if (.not. has_pending) exit worker_source = local_status%MPI_SOURCE if (worker_fragment_map(worker_source) == 0) then call logger%error("Received result from worker "//to_char(worker_source)// & " but no fragment was assigned!") call abort_comm(ctx%resources%mpi_comms%world_comm, 1) end if call result_irecv(results(worker_fragment_map(worker_source)), ctx%resources%mpi_comms%node_comm, worker_source, & TAG_WORKER_SCALAR_RESULT, req) call wait(req) if (results(worker_fragment_map(worker_source))%has_error) then call logger%error("Fragment "//to_char(worker_fragment_map(worker_source))// & " calculation failed: "// & results(worker_fragment_map(worker_source))%error%get_message()) call abort_comm(ctx%resources%mpi_comms%world_comm, 1) end if worker_fragment_map(worker_source) = 0 results_received = results_received + 1 if (mod(results_received, max(1_int64, ctx%total_fragments/10)) == 0 .or. & results_received == ctx%total_fragments) then call logger%info(" Processed "//to_char(results_received)//"/"// & to_char(ctx%total_fragments)//" fragments ["// & to_char(coord_timer%get_elapsed_time())//" s]") end if end do end subroutine handle_local_worker_results