16 #ifndef dealii_mpi_mpi_remote_point_evaluation_h
17 #define dealii_mpi_mpi_remote_point_evaluation_h
33 template <
int dim,
int spacedim>
52 template <
int dim,
int spacedim = dim>
109 DistributedComputePointLocationsInternal<dim, spacedim> &data,
121 std::vector<std::pair<int, int>>
cells;
156 template <
typename T>
159 std::vector<T> &output,
160 std::vector<T> &buffer,
162 &evaluation_function)
const;
168 template <
typename T>
172 &evaluation_function)
const;
182 template <
typename T>
185 const std::vector<T> &input,
186 std::vector<T> &buffer,
188 &evaluation_function)
const;
194 template <
typename T>
197 const std::vector<T> &input,
199 &evaluation_function)
const;
205 const std::vector<unsigned int> &
352 template <
int dim,
int spacedim>
353 template <
typename T>
356 std::vector<T> &output,
357 std::vector<T> &buffer,
359 &evaluation_function)
const
361 #ifndef DEAL_II_WITH_MPI
365 (void)evaluation_function;
370 const unsigned int my_rank =
374 output.resize(point_ptrs.back());
375 buffer.resize(
std::max(send_permutation.size() * 2,
376 point_ptrs.back() + send_permutation.size()));
379 ArrayView<T> buffer_eval(buffer.data(), send_permutation.size());
382 ArrayView<T> buffer_comm(buffer.data() + send_permutation.size(),
383 send_permutation.size());
386 evaluation_function(buffer_eval, cell_data);
392 const auto my_rank_local_recv_ptr =
393 std::find(recv_ranks.begin(), recv_ranks.end(), my_rank);
395 if (my_rank_local_recv_ptr != recv_ranks.end())
398 std::distance(recv_ranks.begin(), my_rank_local_recv_ptr);
399 my_rank_local_send = std::distance(send_ranks.begin(),
400 std::find(send_ranks.begin(),
405 for (
unsigned int i = 0; i < send_permutation.size(); ++i)
407 const unsigned int send_index = send_permutation[i];
411 (send_ptrs[my_rank_local_send] <= send_index &&
412 send_index < send_ptrs[my_rank_local_send + 1]))
413 output[recv_permutation[send_index - send_ptrs[my_rank_local_send] +
414 recv_ptrs[my_rank_local_recv]]] =
417 buffer_comm[send_index] = buffer_eval[i];
421 std::vector<std::vector<char>> send_buffer;
422 send_buffer.reserve(send_ranks.size());
424 std::vector<MPI_Request> send_requests;
425 send_requests.reserve(send_ranks.size());
427 for (
unsigned int i = 0; i < send_ranks.size(); ++i)
429 if (send_ranks[i] == my_rank)
432 send_requests.emplace_back(MPI_Request());
435 std::vector<T>(buffer_comm.
begin() + send_ptrs[i],
436 buffer_comm.
begin() + send_ptrs[i + 1]),
439 const int ierr = MPI_Isend(send_buffer.back().data(),
440 send_buffer.back().size(),
445 &send_requests.back());
450 std::vector<char> buffer_char;
452 for (
unsigned int i = 0; i < recv_ranks.size(); ++i)
454 if (recv_ranks[i] == my_rank)
459 int ierr = MPI_Probe(MPI_ANY_SOURCE,
466 ierr = MPI_Get_count(&status, MPI_CHAR, &message_length);
469 buffer_char.resize(message_length);
471 ierr = MPI_Recv(buffer_char.data(),
482 Utilities::unpack<std::vector<T>>(buffer_char,
false);
486 std::find(recv_ranks.begin(), recv_ranks.end(), status.MPI_SOURCE);
490 const unsigned int j = std::distance(recv_ranks.begin(), ptr);
494 for (
unsigned int i = recv_ptrs[j], c = 0; i < recv_ptrs[j + 1];
496 output[recv_permutation[i]] = buffer[c];
500 if (!send_requests.empty())
502 const int ierr = MPI_Waitall(send_requests.size(),
503 send_requests.data(),
504 MPI_STATUSES_IGNORE);
511 template <
int dim,
int spacedim>
512 template <
typename T>
516 &evaluation_function)
const
518 std::vector<T> output;
519 std::vector<T> buffer;
521 this->evaluate_and_process(output, buffer, evaluation_function);
528 template <
int dim,
int spacedim>
529 template <
typename T>
532 const std::vector<T> &input,
533 std::vector<T> &buffer,
535 &evaluation_function)
const
537 #ifndef DEAL_II_WITH_MPI
541 (void)evaluation_function;
546 const unsigned int my_rank =
550 std::vector<unsigned int> recv_permutation_inv(recv_permutation.size());
551 for (
unsigned int c = 0; c < recv_permutation.size(); ++c)
552 recv_permutation_inv[recv_permutation[c]] = c;
554 std::vector<unsigned int> send_permutation_inv(send_permutation.size());
555 for (
unsigned int c = 0; c < send_permutation.size(); ++c)
556 send_permutation_inv[send_permutation[c]] = c;
559 const auto &point_ptrs = this->get_point_ptrs();
561 buffer.resize(
std::max(send_permutation.size() * 2,
562 point_ptrs.back() + send_permutation.size()));
565 ArrayView<T> buffer_eval(buffer.data(), send_permutation.size());
568 ArrayView<T> buffer_comm(buffer.data() + send_permutation.size(),
575 const auto my_rank_local_recv_ptr =
576 std::find(recv_ranks.begin(), recv_ranks.end(), my_rank);
578 if (my_rank_local_recv_ptr != recv_ranks.end())
581 std::distance(recv_ranks.begin(), my_rank_local_recv_ptr);
582 my_rank_local_send = std::distance(send_ranks.begin(),
583 std::find(send_ranks.begin(),
588 for (
unsigned int i = 0, c = 0; i < point_ptrs.size() - 1; ++i)
590 const auto n_entries = point_ptrs[i + 1] - point_ptrs[i];
592 for (
unsigned int j = 0; j < n_entries; ++j, ++c)
594 const unsigned int recv_index = recv_permutation_inv[c];
598 (recv_ptrs[my_rank_local_recv] <= recv_index &&
599 recv_index < recv_ptrs[my_rank_local_recv + 1]))
600 buffer_eval[send_permutation_inv
601 [recv_index - recv_ptrs[my_rank_local_recv] +
602 send_ptrs[my_rank_local_send]]] = input[i];
604 buffer_comm[recv_index] = input[i];
609 std::vector<std::vector<char>> send_buffer;
610 send_buffer.reserve(recv_ranks.size());
612 std::vector<MPI_Request> send_requests;
613 send_requests.reserve(recv_ranks.size());
615 for (
unsigned int i = 0; i < recv_ranks.size(); ++i)
617 if (recv_ranks[i] == my_rank)
620 send_requests.push_back(MPI_Request());
623 std::vector<T>(buffer_comm.
begin() + recv_ptrs[i],
624 buffer_comm.
begin() + recv_ptrs[i + 1]),
627 const int ierr = MPI_Isend(send_buffer.back().data(),
628 send_buffer.back().size(),
633 &send_requests.back());
638 std::vector<char> recv_buffer;
640 for (
unsigned int i = 0; i < send_ranks.size(); ++i)
642 if (send_ranks[i] == my_rank)
647 int ierr = MPI_Probe(MPI_ANY_SOURCE,
654 ierr = MPI_Get_count(&status, MPI_CHAR, &message_length);
657 recv_buffer.resize(message_length);
659 ierr = MPI_Recv(recv_buffer.data(),
669 const auto recv_buffer_unpacked =
670 Utilities::unpack<std::vector<T>>(recv_buffer,
false);
674 std::find(send_ranks.begin(), send_ranks.end(), status.MPI_SOURCE);
678 const unsigned int j = std::distance(send_ranks.begin(), ptr);
681 send_ptrs[j + 1] - send_ptrs[j]);
683 for (
unsigned int i = send_ptrs[j], c = 0; i < send_ptrs[j + 1];
685 buffer_eval[send_permutation_inv[i]] = recv_buffer_unpacked[c];
688 if (!send_requests.empty())
690 const int ierr = MPI_Waitall(send_requests.size(),
691 send_requests.data(),
692 MPI_STATUSES_IGNORE);
697 evaluation_function(buffer_eval, cell_data);
703 template <
int dim,
int spacedim>
704 template <
typename T>
707 const std::vector<T> &input,
709 &evaluation_function)
const
711 std::vector<T> buffer;
712 this->process_and_evaluate(input, buffer, evaluation_function);
virtual MPI_Comm get_communicator() const
RemotePointEvaluation(const double tolerance=1e-6, const bool enforce_unique_mapping=false, const unsigned int rtree_level=0, const std::function< std::vector< bool >()> &marked_vertices={})
const std::function< std::vector< bool >)> marked_vertices
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
bool all_points_found() const
const bool enforce_unique_mapping
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 unsigned int rtree_level
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
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
void reinit(const std::vector< Point< spacedim >> &points, const Triangulation< dim, spacedim > &tria, const Mapping< dim, spacedim > &mapping)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNeedsMPI()
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
#define AssertThrowMPI(error_code)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
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
std::vector< std::pair< int, int > > cells
std::vector< unsigned int > reference_point_ptrs
std::vector< Point< dim > > reference_point_values
const ::Triangulation< dim, spacedim > & tria