15#ifndef dealii_mpi_mpi_remote_point_evaluation_h
16#define dealii_mpi_mpi_remote_point_evaluation_h
32 template <
int dim,
int spacedim>
51 template <
int dim,
int spacedim = dim>
125 DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
126 "Use the constructor with AdditionalData struct.")
128 const
double tolerance,
129 const
bool enforce_unique_mapping = false,
130 const
unsigned int rtree_level = 0,
131 const
std::function<
std::vector<
bool>()> &marked_vertices = {});
165 DistributedComputePointLocationsInternal<dim, spacedim> &data,
226 template <
typename T>
234 std::vector<std::pair<int, int>>
cells;
276 template <
typename T>
279 std::vector<T> &output,
280 std::vector<T> &buffer,
282 &evaluation_function)
const;
288 template <
typename T>
292 &evaluation_function)
const;
302 template <
typename T>
305 const std::vector<T> &input,
306 std::vector<T> &buffer,
308 &evaluation_function)
const;
314 template <
typename T>
317 const std::vector<T> &input,
319 &evaluation_function)
const;
325 const std::vector<unsigned int> &
455 template <
int dim,
int spacedim>
456 template <
typename T>
459 const unsigned int cell,
469 template <
int dim,
int spacedim>
470 template <
typename T>
473 std::vector<T> &output,
474 std::vector<T> &buffer,
476 &evaluation_function)
const
478#ifndef DEAL_II_WITH_MPI
482 (void)evaluation_function;
487 const unsigned int my_rank =
503 evaluation_function(buffer_eval, *
cell_data);
509 const auto my_rank_local_recv_ptr =
512 if (my_rank_local_recv_ptr !=
recv_ranks.end())
515 std::distance(
recv_ranks.begin(), my_rank_local_recv_ptr);
516 my_rank_local_send = std::distance(
send_ranks.begin(),
528 (
send_ptrs[my_rank_local_send] <= send_index &&
529 send_index <
send_ptrs[my_rank_local_send + 1]))
534 buffer_comm[send_index] = buffer_eval[i];
538 std::vector<std::vector<char>> send_buffer;
541 std::vector<MPI_Request> send_requests;
544 for (
unsigned int i = 0; i <
send_ranks.size(); ++i)
549 send_requests.emplace_back(MPI_Request());
556 const int ierr = MPI_Isend(send_buffer.back().data(),
557 send_buffer.back().size(),
562 &send_requests.back());
567 std::vector<char> buffer_char;
569 for (
unsigned int i = 0; i <
recv_ranks.size(); ++i)
576 int ierr = MPI_Probe(MPI_ANY_SOURCE,
583 ierr = MPI_Get_count(&status, MPI_CHAR, &message_length);
586 buffer_char.resize(message_length);
588 ierr = MPI_Recv(buffer_char.data(),
599 Utilities::unpack<std::vector<T>>(buffer_char,
false);
607 const unsigned int j = std::distance(
recv_ranks.begin(), ptr);
617 if (!send_requests.empty())
619 const int ierr = MPI_Waitall(send_requests.size(),
620 send_requests.data(),
621 MPI_STATUSES_IGNORE);
628 template <
int dim,
int spacedim>
629 template <
typename T>
633 &evaluation_function)
const
635 std::vector<T> output;
636 std::vector<T> buffer;
645 template <
int dim,
int spacedim>
646 template <
typename T>
649 const std::vector<T> &input,
650 std::vector<T> &buffer,
652 &evaluation_function)
const
654#ifndef DEAL_II_WITH_MPI
658 (void)evaluation_function;
663 const unsigned int my_rank =
692 const auto my_rank_local_recv_ptr =
695 if (my_rank_local_recv_ptr !=
recv_ranks.end())
698 std::distance(
recv_ranks.begin(), my_rank_local_recv_ptr);
699 my_rank_local_send = std::distance(
send_ranks.begin(),
705 for (
unsigned int i = 0, c = 0; i <
point_ptrs.size() - 1; ++i)
709 for (
unsigned int j = 0; j < n_entries; ++j, ++c)
711 const unsigned int recv_index = recv_permutation_inv[c];
715 (
recv_ptrs[my_rank_local_recv] <= recv_index &&
716 recv_index <
recv_ptrs[my_rank_local_recv + 1]))
717 buffer_eval[send_permutation_inv
718 [recv_index -
recv_ptrs[my_rank_local_recv] +
719 send_ptrs[my_rank_local_send]]] = input[i];
721 buffer_comm[recv_index] = input[i];
726 std::vector<std::vector<char>> send_buffer;
729 std::vector<MPI_Request> send_requests;
732 for (
unsigned int i = 0; i <
recv_ranks.size(); ++i)
737 send_requests.push_back(MPI_Request());
744 const int ierr = MPI_Isend(send_buffer.back().data(),
745 send_buffer.back().size(),
750 &send_requests.back());
755 std::vector<char> recv_buffer;
757 for (
unsigned int i = 0; i <
send_ranks.size(); ++i)
764 int ierr = MPI_Probe(MPI_ANY_SOURCE,
771 ierr = MPI_Get_count(&status, MPI_CHAR, &message_length);
774 recv_buffer.resize(message_length);
776 ierr = MPI_Recv(recv_buffer.data(),
786 const auto recv_buffer_unpacked =
787 Utilities::unpack<std::vector<T>>(recv_buffer,
false);
795 const unsigned int j = std::distance(
send_ranks.begin(), ptr);
802 buffer_eval[send_permutation_inv[i]] = recv_buffer_unpacked[c];
805 if (!send_requests.empty())
807 const int ierr = MPI_Waitall(send_requests.size(),
808 send_requests.data(),
809 MPI_STATUSES_IGNORE);
814 evaluation_function(buffer_eval, *
cell_data);
820 template <
int dim,
int spacedim>
821 template <
typename T>
824 const std::vector<T> &input,
826 &evaluation_function)
const
828 std::vector<T> buffer;
Abstract base class for mapping classes.
virtual MPI_Comm get_communicator() const
ArrayView< T > get_data_view(const unsigned int cell, const ArrayView< T > &values) const
std::vector< std::pair< int, int > > cells
std_cxx20::ranges::iota_view< unsigned int, unsigned int > cell_indices() const
ArrayView< const Point< dim > > get_unit_points(const unsigned int cell) const
Triangulation< dim, spacedim >::active_cell_iterator get_active_cell_iterator(const unsigned int cell) const
const Triangulation< dim, spacedim > & triangulation
std::vector< unsigned int > reference_point_ptrs
std::vector< Point< dim > > reference_point_values
void evaluate_and_process(std::vector< T > &output, std::vector< T > &buffer, const std::function< void(const ArrayView< T > &, const CellData &)> &evaluation_function) const
std::vector< unsigned int > send_ptrs
std::unique_ptr< CellData > cell_data
bool all_points_found() const
boost::signals2::connection tria_signal
std::vector< unsigned int > point_ptrs
std::vector< unsigned int > recv_permutation
void process_and_evaluate(const std::vector< T > &input, const std::function< void(const ArrayView< const T > &, const CellData &)> &evaluation_function) const
const Triangulation< dim, spacedim > & get_triangulation() const
SmartPointer< const Mapping< dim, spacedim > > mapping
void process_and_evaluate(const std::vector< T > &input, std::vector< T > &buffer, const std::function< void(const ArrayView< const T > &, const CellData &)> &evaluation_function) const
const std::vector< unsigned int > & get_point_ptrs() const
const AdditionalData additional_data
bool is_map_unique() const
std::vector< unsigned int > recv_ranks
std::vector< unsigned int > recv_ptrs
const Mapping< dim, spacedim > & get_mapping() const
const CellData & get_cell_data() const
std::vector< unsigned int > send_ranks
bool all_points_found_flag
void reinit(const std::vector< Point< spacedim > > &points, const Triangulation< dim, spacedim > &tria, const Mapping< dim, spacedim > &mapping)
SmartPointer< const Triangulation< dim, spacedim > > tria
std::vector< T > evaluate_and_process(const std::function< void(const ArrayView< T > &, const CellData &)> &evaluation_function) const
bool point_found(const unsigned int i) const
std::vector< unsigned int > send_permutation
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcNeedsMPI()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertThrowMPI(error_code)
#define AssertIndexRange(index, range)
unsigned int this_mpi_process(const MPI_Comm mpi_communicator)
size_t pack(const T &object, std::vector< char > &dest_buffer, const bool allow_compression=true)
static const unsigned int invalid_unsigned_int
boost::integer_range< IncrementableType > iota_view
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
std::function< std::vector< bool >()> marked_vertices
bool enforce_unique_mapping