fill send buffers and start sending
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(tem_communication_type), | intent(inout) | :: | send | |||
type(tem_communication_type), | intent(inout) | :: | recv | |||
real(kind=rk), | intent(inout) | :: | state(*) | |||
integer, | intent(in) | :: | message_flag | |||
real(kind=rk), | intent(in), | optional | :: | send_state(*) | ||
integer, | intent(in) | :: | comm |
mpi communicator |
subroutine comm_isend_irecv_overlap_real( send, recv, state, & & message_flag, send_state, comm ) ! -------------------------------------------------------------------- ! type( tem_communication_type), intent(inout) :: send, recv real(kind=rk), intent(inout) :: state(*) integer, intent(in) :: message_flag real(kind=rk), intent(in), optional :: send_state(*) !< data to send !> mpi communicator integer, intent(in) :: comm ! -------------------------------------------------------------------- ! integer :: iproc, ival ! counter for neigbor mees integer :: finproc ! integer :: rq_handle( recv%nprocs + send%nprocs ) integer :: status( mpi_status_size, max(send%nprocs,1) ) integer :: ierr ! error flag integer :: nsendvals, nrecvvals ! -------------------------------------------------------------------- ! ! pre-post irecv do iproc = 1, recv%nprocs call mpi_irecv( & & recv%buf_real( iproc )%val, & ! buffer & recv%buf_real( iproc )%nvals, & ! count & rk_mpi, & ! data type & recv%proc(iproc), & ! target me & message_flag, & ! tag & comm, & ! communicator & recv%rqhandle(iproc), & ! handle & ierr ) ! error status end do !> fill send buffers and start sending do iproc = 1, send%nprocs nsendvals = send%buf_real( iproc )%nvals if (present(send_state)) then !$nec ivdep do ival = 1, nsendvals send%buf_real( iproc )%val( ival ) & & = send_state( send%buf_real( iproc )%pos( ival ) ) end do else !$nec ivdep do ival = 1, nsendvals send%buf_real( iproc )%val( ival ) & & = state( send%buf_real( iproc )%pos( ival ) ) end do end if call mpi_isend( & & send%buf_real( iproc )%val, & ! buffer & send%buf_real( iproc )%nvals, & ! count & rk_mpi, & ! data type & send%proc(iproc), & ! target & message_flag, & ! tag & comm, & ! comm & send%rqhandle( iproc ), & ! handle & ierr ) ! error status end do do iproc = 1, recv%nprocs ! wait for any of receive buffer to be ready call mpi_waitany( & & recv%nprocs, & ! count & recv%rqhandle, & ! request handles & finproc, & ! process that finished & status(:,1), & ! mpi status & ierr ) ! error status nrecvvals = recv%buf_real( finproc )%nvals !$nec ivdep do ival = 1, nrecvvals ! write the values from the recv me to the state array ! to positions given in recvme%pos state( recv%buf_real( finproc )%pos( ival ) ) & & = recv%buf_real( finproc )%val( ival ) end do end do if (send%nprocs > 0) then ! wait for send buffer to be ready call mpi_waitall( & & send%nprocs, & ! total number of comm.'s to wait for & send%rqhandle, & ! request handles & status, & ! mpi status & ierr ) ! error status end if end subroutine comm_isend_irecv_overlap_real