flow_exchange_cell_matrix Subroutine

public subroutine flow_exchange_cell_matrix(flow, field)

Uses

  • proc~~flow_exchange_cell_matrix~~UsesGraph proc~flow_exchange_cell_matrix mod_mpi_flow::flow_exchange_cell_matrix module~mod_profiler mod_profiler proc~flow_exchange_cell_matrix->module~mod_profiler iso_fortran_env iso_fortran_env module~mod_profiler->iso_fortran_env module~mod_kinds mod_kinds module~mod_profiler->module~mod_kinds mpi_f08 mpi_f08 module~mod_profiler->mpi_f08 module~mod_kinds->iso_fortran_env

Exchanges owned cell matrix values to ranks that keep them as ghosts.

Arguments

Type IntentOptional Attributes Name
type(flow_mpi_t), intent(inout) :: flow
real(kind=rk), intent(inout) :: field(:,:)

Calls

proc~~flow_exchange_cell_matrix~~CallsGraph proc~flow_exchange_cell_matrix mod_mpi_flow::flow_exchange_cell_matrix mpi_irecv mpi_irecv proc~flow_exchange_cell_matrix->mpi_irecv mpi_isend mpi_isend proc~flow_exchange_cell_matrix->mpi_isend mpi_waitall mpi_waitall proc~flow_exchange_cell_matrix->mpi_waitall proc~check_mpi~2 mod_mpi_flow::check_mpi proc~flow_exchange_cell_matrix->proc~check_mpi~2 proc~fatal_error mod_kinds::fatal_error proc~flow_exchange_cell_matrix->proc~fatal_error proc~profiler_start mod_profiler::profiler_start proc~flow_exchange_cell_matrix->proc~profiler_start proc~profiler_stop mod_profiler::profiler_stop proc~flow_exchange_cell_matrix->proc~profiler_stop proc~check_mpi~2->proc~fatal_error mpi_wtime mpi_wtime proc~profiler_start->mpi_wtime proc~find_or_create_timer mod_profiler::find_or_create_timer proc~profiler_start->proc~find_or_create_timer proc~profiler_stop->mpi_wtime proc~profiler_stop->proc~find_or_create_timer proc~record_edge mod_profiler::record_edge proc~profiler_stop->proc~record_edge

Called by

proc~~flow_exchange_cell_matrix~~CalledByGraph proc~flow_exchange_cell_matrix mod_mpi_flow::flow_exchange_cell_matrix proc~advance_projection_step mod_flow_projection::advance_projection_step proc~advance_projection_step->proc~flow_exchange_cell_matrix proc~advance_species_transport mod_species::advance_species_transport proc~advance_species_transport->proc~flow_exchange_cell_matrix proc~update_transport_properties mod_transport_properties::update_transport_properties proc~update_transport_properties->proc~flow_exchange_cell_matrix program~lowmach_react_hex lowmach_react_hex program~lowmach_react_hex->proc~advance_projection_step program~lowmach_react_hex->proc~advance_species_transport program~lowmach_react_hex->proc~update_transport_properties

Source Code

   subroutine flow_exchange_cell_matrix(flow, field)
      use mod_profiler, only : profiler_start, profiler_stop
      type(flow_mpi_t), intent(inout) :: flow
      real(rk), intent(inout) :: field(:,:)
      integer, parameter :: cell_matrix_halo_tag = 9282
      integer :: ncomp, i, j, k, nreq, ierr, offset, count, pos

      if (.not. allocated(flow%cell_sendbuf)) return
      if (flow%ncell_recv_ranks + flow%ncell_send_ranks == 0) return

      ncomp = size(field, 1)
      if (ncomp > flow%cell_halo_max_components) then
         call fatal_error('mpi_flow', 'cell matrix halo component count exceeds cached buffer size')
      end if

      do i = 1, flow%ncell_send_ranks
         offset = flow%cell_send_displs(i)
         count = flow%cell_send_counts(i)
         do j = 1, count
            pos = (offset + j - 1) * ncomp
            do k = 1, ncomp
               flow%cell_sendbuf(pos + k) = field(k, flow%cell_send_cells(offset + j))
            end do
         end do
      end do

      call profiler_start('MPI_Communication')
      nreq = 0
      do i = 1, flow%ncell_recv_ranks
         offset = flow%cell_recv_displs(i) * ncomp
         count = flow%cell_recv_counts(i) * ncomp
         nreq = nreq + 1
         call MPI_Irecv(flow%cell_recvbuf(offset + 1), count, MPI_DOUBLE_PRECISION, &
                        flow%cell_recv_ranks(i), cell_matrix_halo_tag, flow%comm, flow%cell_requests(nreq), ierr)
         call check_mpi(ierr, 'cell matrix halo irecv')
      end do
      do i = 1, flow%ncell_send_ranks
         offset = flow%cell_send_displs(i) * ncomp
         count = flow%cell_send_counts(i) * ncomp
         nreq = nreq + 1
         call MPI_Isend(flow%cell_sendbuf(offset + 1), count, MPI_DOUBLE_PRECISION, &
                        flow%cell_send_ranks(i), cell_matrix_halo_tag, flow%comm, flow%cell_requests(nreq), ierr)
         call check_mpi(ierr, 'cell matrix halo isend')
      end do
      call MPI_Waitall(nreq, flow%cell_requests(1:nreq), MPI_STATUSES_IGNORE, ierr)
      call check_mpi(ierr, 'cell matrix halo waitall')
      call profiler_stop('MPI_Communication')

      do i = 1, size(flow%cell_recv_cells)
         pos = (i - 1) * ncomp
         do k = 1, ncomp
            field(k, flow%cell_recv_cells(i)) = flow%cell_recvbuf(pos + k)
         end do
      end do
   end subroutine flow_exchange_cell_matrix