16 #include <deal.II/base/std_cxx14/memory.h> 18 #include <deal.II/grid/grid_tools.h> 19 #include <deal.II/grid/grid_tools_cache.h> 21 #include <deal.II/particles/particle_handler.h> 25 DEAL_II_NAMESPACE_OPEN
27 #ifdef DEAL_II_WITH_P4EST 31 template <
int dim,
int spacedim>
36 , global_number_of_particles(0)
37 , global_max_particles_per_cell(0)
38 , next_free_particle_index(0)
43 , handle(
numbers::invalid_unsigned_int)
48 template <
int dim,
int spacedim>
52 const unsigned int n_properties)
53 : triangulation(&triangulation, typeid(*this).name())
54 , mapping(&mapping, typeid(*this).name())
57 , global_number_of_particles(0)
58 , global_max_particles_per_cell(0)
59 , next_free_particle_index(0)
64 , handle(
numbers::invalid_unsigned_int)
69 template <
int dim,
int spacedim>
75 const unsigned int n_properties)
77 triangulation = &new_triangulation;
78 mapping = &new_mapping;
81 property_pool = std_cxx14::make_unique<PropertyPool>(n_properties);
86 template <
int dim,
int spacedim>
91 global_number_of_particles = 0;
92 next_free_particle_index = 0;
93 global_max_particles_per_cell = 0;
98 template <
int dim,
int spacedim>
107 template <
int dim,
int spacedim>
112 unsigned int local_max_particles_per_cell = 0;
113 unsigned int current_particles_per_cell = 0;
115 triangulation->begin_active();
119 locally_highest_index =
120 std::max(locally_highest_index, particle->get_id());
122 if (particle->get_surrounding_cell(*triangulation) != current_cell)
124 current_particles_per_cell = 0;
125 current_cell = particle->get_surrounding_cell(*triangulation);
128 ++current_particles_per_cell;
129 local_max_particles_per_cell =
130 std::max(local_max_particles_per_cell, current_particles_per_cell);
133 global_number_of_particles =
134 ::Utilities::MPI::sum(particles.size(),
135 triangulation->get_communicator());
136 next_free_particle_index =
137 ::Utilities::MPI::max(locally_highest_index,
138 triangulation->get_communicator()) +
140 global_max_particles_per_cell =
141 ::Utilities::MPI::max(local_max_particles_per_cell,
142 triangulation->get_communicator());
147 template <
int dim,
int spacedim>
156 template <
int dim,
int spacedim>
165 template <
int dim,
int spacedim>
174 template <
int dim,
int spacedim>
183 template <
int dim,
int spacedim>
192 template <
int dim,
int spacedim>
201 template <
int dim,
int spacedim>
210 template <
int dim,
int spacedim>
219 template <
int dim,
int spacedim>
226 ->particles_in_cell(cell);
231 template <
int dim,
int spacedim>
236 const internal::LevelInd level_index =
237 std::make_pair<int, int>(cell->level(), cell->index());
239 if (cell->is_ghost())
241 const auto particles_in_cell = ghost_particles.equal_range(level_index);
242 return boost::make_iterator_range(
247 const auto particles_in_cell = particles.equal_range(level_index);
248 return boost::make_iterator_range(
255 template <
int dim,
int spacedim>
260 particles.erase(particle->particle);
265 template <
int dim,
int spacedim>
271 typename std::multimap<internal::LevelInd,
274 std::make_pair(internal::LevelInd(cell->level(), cell->index()),
278 particle_it->set_property_pool(*property_pool);
281 for (
unsigned int n = 0; n < particle.
get_properties().size(); ++n)
289 template <
int dim,
int spacedim>
296 for (
auto particle = new_particles.begin(); particle != new_particles.end();
300 std::make_pair(internal::LevelInd(particle->first->level(),
301 particle->first->index()),
304 update_cached_numbers();
309 template <
int dim,
int spacedim>
314 update_cached_numbers();
321 get_next_free_particle_index();
324 # ifdef DEAL_II_WITH_MPI 326 const int ierr = MPI_Scan(&particles_to_add_locally,
329 DEAL_II_PARTICLE_INDEX_MPI_TYPE,
331 triangulation->get_communicator());
333 local_start_index -= particles_to_add_locally;
336 local_start_index += local_next_particle_index;
341 auto &cells = std::get<0>(point_locations);
342 auto &local_positions = std::get<1>(point_locations);
343 auto &index_map = std::get<2>(point_locations);
345 if (cells.size() == 0)
349 particles.find(std::make_pair(cells[0]->level(), cells[0]->index()));
350 for (
unsigned int i = 0; i < cells.size(); ++i)
352 internal::LevelInd current_cell(cells[i]->level(), cells[i]->index());
353 for (
unsigned int p = 0; p < local_positions[i].size(); ++p)
355 hint = particles.insert(
357 std::make_pair(current_cell,
359 local_positions[i][p],
365 update_cached_numbers();
370 template <
int dim,
int spacedim>
374 return global_number_of_particles;
379 template <
int dim,
int spacedim>
383 return global_max_particles_per_cell;
388 template <
int dim,
int spacedim>
392 return particles.size();
397 template <
int dim,
int spacedim>
401 return property_pool->n_properties_per_slot();
406 template <
int dim,
int spacedim>
412 const internal::LevelInd found_cell =
413 std::make_pair<int, int>(cell->level(), cell->index());
415 if (cell->is_locally_owned())
416 return particles.count(found_cell);
417 else if (cell->is_ghost())
418 return ghost_particles.count(found_cell);
419 else if (cell->is_artificial())
427 template <
int dim,
int spacedim>
431 return next_free_particle_index;
436 template <
int dim,
int spacedim>
440 return *property_pool;
458 compare_particle_association(
459 const unsigned int a,
460 const unsigned int b,
464 const double scalar_product_a = center_directions[a] * particle_direction;
465 const double scalar_product_b = center_directions[b] * particle_direction;
470 return (scalar_product_a > scalar_product_b);
476 template <
int dim,
int spacedim>
488 std::vector<particle_iterator> particles_out_of_cell;
489 particles_out_of_cell.reserve(n_locally_owned_particles());
495 it->get_surrounding_cell(*triangulation);
500 mapping->transform_real_to_unit_cell(cell, it->get_location());
503 it->set_reference_location(p_unit);
508 particles_out_of_cell.push_back(it);
514 particles_out_of_cell.push_back(it);
524 std::vector<std::pair<internal::LevelInd, Particle<dim, spacedim>>>
526 std::map<types::subdomain_id, std::vector<particle_iterator>>
530 std::vector<typename Triangulation<dim, spacedim>::active_cell_iterator>>
538 using vector_size =
typename std::vector<particle_iterator>::size_type;
539 sorted_particles.reserve(
540 static_cast<vector_size>(particles_out_of_cell.size() * 1.25));
542 const std::set<types::subdomain_id> ghost_owners =
543 triangulation->ghost_owners();
545 for (
const auto ghost_owner : ghost_owners)
546 moved_particles[ghost_owner].reserve(
547 static_cast<vector_size>(particles_out_of_cell.size() * 0.25));
548 for (
const auto ghost_owner : ghost_owners)
549 moved_cells[ghost_owner].reserve(
550 static_cast<vector_size>(particles_out_of_cell.size() * 0.25));
555 std::set<typename Triangulation<dim, spacedim>::active_cell_iterator>>
559 const std::vector<std::vector<Tensor<1, spacedim>>>
560 vertex_to_cell_centers(
564 std::vector<unsigned int> neighbor_permutation;
567 typename std::vector<particle_iterator>::iterator
568 it = particles_out_of_cell.begin(),
569 end_particle = particles_out_of_cell.end();
571 for (; it != end_particle; ++it)
575 bool found_cell =
false;
580 current_cell = (*it)->get_surrounding_cell(*triangulation);
582 const unsigned int closest_vertex =
583 GridTools::find_closest_vertex_of_cell<dim, spacedim>(
584 current_cell, (*it)->get_location());
586 (*it)->get_location() - current_cell->vertex(closest_vertex);
587 vertex_to_particle /= vertex_to_particle.
norm();
589 const unsigned int closest_vertex_index =
590 current_cell->vertex_index(closest_vertex);
591 const unsigned int n_neighbor_cells =
592 vertex_to_cells[closest_vertex_index].size();
594 neighbor_permutation.resize(n_neighbor_cells);
595 for (
unsigned int i = 0; i < n_neighbor_cells; ++i)
596 neighbor_permutation[i] = i;
598 std::sort(neighbor_permutation.begin(),
599 neighbor_permutation.end(),
600 std::bind(&compare_particle_association<spacedim>,
601 std::placeholders::_1,
602 std::placeholders::_2,
603 std::cref(vertex_to_particle),
605 vertex_to_cell_centers[closest_vertex_index])));
609 for (
unsigned int i = 0; i < n_neighbor_cells; ++i)
613 typename std::set<typename Triangulation<dim, spacedim>::
614 active_cell_iterator>::const_iterator
615 cell = vertex_to_cells[closest_vertex_index].begin();
617 std::advance(cell, neighbor_permutation[i]);
619 mapping->transform_real_to_unit_cell(*cell,
620 (*it)->get_location());
623 current_cell = *cell;
624 current_reference_position = p_unit;
640 const std::pair<const typename Triangulation<dim, spacedim>::
641 active_cell_iterator,
643 current_cell_and_position =
644 GridTools::find_active_cell_around_point<>(
645 *mapping, *triangulation, (*it)->get_location());
646 current_cell = current_cell_and_position.first;
647 current_reference_position = current_cell_and_position.second;
649 catch (GridTools::ExcPointNotFound<spacedim> &)
659 (*it)->set_reference_location(current_reference_position);
663 if (current_cell->is_locally_owned())
665 sorted_particles.push_back(
666 std::make_pair(internal::LevelInd(current_cell->level(),
667 current_cell->index()),
668 (*it)->particle->second));
672 moved_particles[current_cell->subdomain_id()].push_back(*it);
673 moved_cells[current_cell->subdomain_id()].push_back(current_cell);
680 std::multimap<internal::LevelInd, Particle<dim, spacedim>>
681 sorted_particles_map;
684 # ifdef DEAL_II_WITH_MPI 686 triangulation->get_communicator()) > 1)
687 send_recv_particles(moved_particles, sorted_particles_map, moved_cells);
690 sorted_particles_map.insert(sorted_particles.begin(),
691 sorted_particles.end());
693 for (
unsigned int i = 0; i < particles_out_of_cell.size(); ++i)
694 remove_particle(particles_out_of_cell[i]);
696 particles.insert(sorted_particles_map.begin(), sorted_particles_map.end());
697 update_cached_numbers();
702 template <
int dim,
int spacedim>
708 triangulation->get_communicator()) == 1)
711 # ifdef DEAL_II_WITH_MPI 713 ghost_particles.clear();
715 std::map<types::subdomain_id, std::vector<particle_iterator>>
716 ghost_particles_by_domain;
718 const std::set<types::subdomain_id> ghost_owners =
719 triangulation->ghost_owners();
720 for (
const auto ghost_owner : ghost_owners)
721 ghost_particles_by_domain[ghost_owner].reserve(
722 static_cast<typename std::vector<particle_iterator>::size_type
>(
723 particles.size() * 0.25));
725 std::vector<std::set<unsigned int>> vertex_to_neighbor_subdomain(
726 triangulation->n_vertices());
728 for (
const auto &cell : triangulation->active_cell_iterators())
730 if (cell->is_ghost())
731 for (
unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell;
733 vertex_to_neighbor_subdomain[cell->vertex_index(v)].insert(
734 cell->subdomain_id());
737 for (
const auto &cell : triangulation->active_cell_iterators())
739 if (!cell->is_ghost())
741 std::set<unsigned int> cell_to_neighbor_subdomain;
742 for (
unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell;
745 cell_to_neighbor_subdomain.insert(
746 vertex_to_neighbor_subdomain[cell->vertex_index(v)].begin(),
747 vertex_to_neighbor_subdomain[cell->vertex_index(v)].end());
750 if (cell_to_neighbor_subdomain.size() > 0)
753 particles_in_cell(cell);
755 for (
const auto domain : cell_to_neighbor_subdomain)
757 for (
typename particle_iterator_range::iterator particle =
758 particle_range.begin();
759 particle != particle_range.end();
761 ghost_particles_by_domain[domain].push_back(particle);
767 send_recv_particles(ghost_particles_by_domain, ghost_particles);
773 # ifdef DEAL_II_WITH_MPI 774 template <
int dim,
int spacedim>
787 const std::set<types::subdomain_id> ghost_owners =
788 triangulation->ghost_owners();
789 const std::vector<types::subdomain_id> neighbors(ghost_owners.begin(),
791 const unsigned int n_neighbors = neighbors.size();
793 if (send_cells.size() != 0)
799 particles_to_send.end(),
804 for (
auto send_particles = particles_to_send.begin();
805 send_particles != particles_to_send.end();
807 Assert(ghost_owners.find(send_particles->first) != ghost_owners.end(),
810 unsigned int n_send_particles = 0;
811 for (
auto send_particles = particles_to_send.begin();
812 send_particles != particles_to_send.end();
814 n_send_particles += send_particles->second.size();
820 std::vector<unsigned int> n_send_data(n_neighbors, 0);
821 std::vector<unsigned int> send_offsets(n_neighbors, 0);
822 std::vector<char> send_data;
827 if (n_send_particles > 0)
830 const unsigned int particle_size =
831 begin()->serialized_size_in_bytes() + cellid_size +
832 (size_callback ? size_callback() : 0);
833 send_data.resize(n_send_particles * particle_size);
834 void *data =
static_cast<void *
>(&send_data.front());
837 for (
unsigned int i = 0; i < n_neighbors; ++i)
839 send_offsets[i] =
reinterpret_cast<std::size_t
>(data) -
840 reinterpret_cast<std::size_t>(&send_data.front());
842 for (
unsigned int j = 0;
843 j < particles_to_send.at(neighbors[i]).size();
849 if (send_cells.size() == 0)
851 particles_to_send.at(neighbors[i])[j]->get_surrounding_cell(
854 cell = send_cells.at(neighbors[i])[j];
857 cell->id().template to_binary<dim>();
858 memcpy(data, &cellid, cellid_size);
859 data =
static_cast<char *
>(data) + cellid_size;
861 particles_to_send.at(neighbors[i])[j]->write_data(data);
864 store_callback(particles_to_send.at(neighbors[i])[j], data);
866 n_send_data[i] =
reinterpret_cast<std::size_t
>(data) -
868 reinterpret_cast<std::size_t>(&send_data.front());
873 std::vector<unsigned int> n_recv_data(n_neighbors);
874 std::vector<unsigned int> recv_offsets(n_neighbors);
878 std::vector<MPI_Request> n_requests(2 * n_neighbors);
879 for (
unsigned int i = 0; i < n_neighbors; ++i)
881 const int ierr = MPI_Irecv(&(n_recv_data[i]),
886 triangulation->get_communicator(),
887 &(n_requests[2 * i]));
890 for (
unsigned int i = 0; i < n_neighbors; ++i)
892 const int ierr = MPI_Isend(&(n_send_data[i]),
897 triangulation->get_communicator(),
898 &(n_requests[2 * i + 1]));
902 MPI_Waitall(2 * n_neighbors, n_requests.data(), MPI_STATUSES_IGNORE);
907 unsigned int total_recv_data = 0;
908 for (
unsigned int neighbor_id = 0; neighbor_id < n_neighbors; ++neighbor_id)
910 recv_offsets[neighbor_id] = total_recv_data;
911 total_recv_data += n_recv_data[neighbor_id];
915 std::vector<char> recv_data(total_recv_data);
919 std::vector<MPI_Request> requests(2 * n_neighbors);
920 unsigned int send_ops = 0;
921 unsigned int recv_ops = 0;
923 for (
unsigned int i = 0; i < n_neighbors; ++i)
924 if (n_recv_data[i] > 0)
926 const int ierr = MPI_Irecv(&(recv_data[recv_offsets[i]]),
931 triangulation->get_communicator(),
932 &(requests[send_ops]));
937 for (
unsigned int i = 0; i < n_neighbors; ++i)
938 if (n_send_data[i] > 0)
940 const int ierr = MPI_Isend(&(send_data[send_offsets[i]]),
945 triangulation->get_communicator(),
946 &(requests[send_ops + recv_ops]));
951 MPI_Waitall(send_ops + recv_ops, requests.data(), MPI_STATUSES_IGNORE);
957 const void *recv_data_it =
static_cast<const void *
>(recv_data.data());
959 while (reinterpret_cast<std::size_t>(recv_data_it) -
960 reinterpret_cast<std::size_t
>(recv_data.data()) <
964 memcpy(&binary_cellid, recv_data_it, cellid_size);
965 const CellId id(binary_cellid);
966 recv_data_it =
static_cast<const char *
>(recv_data_it) + cellid_size;
969 id.to_cell(*triangulation);
971 typename std::multimap<internal::LevelInd,
973 recv_particle = received_particles.insert(std::make_pair(
974 internal::LevelInd(cell->level(), cell->index()),
983 AssertThrow(recv_data_it == recv_data.data() + recv_data.size(),
985 "The amount of data that was read into new particles " 986 "does not match the amount of data sent around."));
992 template <
int dim,
int spacedim>
995 const std::function<std::size_t()> & size_callb,
1000 size_callback = size_callb;
1001 store_callback = store_callb;
1002 load_callback = load_callb;
1007 template <
int dim,
int spacedim>
1012 *non_const_triangulation =
1018 update_cached_numbers();
1020 if (global_max_particles_per_cell > 0)
1022 const std::function<std::vector<char>(
1028 std::placeholders::_1,
1029 std::placeholders::_2);
1032 callback_function,
true);
1038 template <
int dim,
int spacedim>
1041 const bool serialization)
1048 *non_const_triangulation =
1056 if (serialization && (global_max_particles_per_cell > 0))
1058 const std::function<std::vector<char>(
1064 std::placeholders::_1,
1065 std::placeholders::_2);
1068 callback_function,
true);
1074 const std::function<void(
1077 const boost::iterator_range<std::vector<char>::const_iterator> &)>
1081 std::placeholders::_1,
1082 std::placeholders::_2,
1083 std::placeholders::_3);
1091 update_cached_numbers();
1097 template <
int dim,
int spacedim>
1103 std::vector<Particle<dim, spacedim>> stored_particles_on_cell;
1112 unsigned int n_particles = 0;
1114 const internal::LevelInd level_index = {cell->level(),
1116 const auto particles_in_cell =
1117 (cell->is_ghost() ? ghost_particles.equal_range(level_index) :
1118 particles.equal_range(level_index));
1120 n_particles = n_particles_in_cell(cell);
1121 stored_particles_on_cell.reserve(n_particles);
1124 particles_in_cell.first,
1125 particles_in_cell.second,
1126 [&stored_particles_on_cell](
1129 stored_particles_on_cell.push_back(particle.second);
1140 unsigned int n_particles = 0;
1142 for (
unsigned int child_index = 0;
1143 child_index < GeometryInfo<dim>::max_children_per_cell;
1147 child = cell->child(child_index);
1148 n_particles += n_particles_in_cell(child);
1151 stored_particles_on_cell.reserve(n_particles);
1153 for (
unsigned int child_index = 0;
1154 child_index < GeometryInfo<dim>::max_children_per_cell;
1158 child = cell->child(child_index);
1159 const internal::LevelInd level_index = {child->level(),
1161 const auto particles_in_cell =
1162 (child->is_ghost() ?
1163 ghost_particles.equal_range(level_index) :
1164 particles.equal_range(level_index));
1167 particles_in_cell.first,
1168 particles_in_cell.second,
1169 [&stored_particles_on_cell](
1172 stored_particles_on_cell.push_back(particle.second);
1189 template <
int dim,
int spacedim>
1194 const boost::iterator_range<std::vector<char>::const_iterator> &data_range)
1198 std::vector<Particle<dim, spacedim>> loaded_particles_on_cell =
1199 Utilities::unpack<std::vector<Particle<dim, spacedim>>>(
1207 for (
auto &particle : loaded_particles_on_cell)
1208 particle.set_property_pool(*property_pool);
1214 auto position_hint = particles.
end();
1215 for (
const auto &particle : loaded_particles_on_cell)
1221 # ifdef DEAL_II_WITH_CXX14 1223 particles.emplace_hint(position_hint,
1224 std::make_pair(cell->level(),
1226 std::move(particle));
1229 particles.insert(position_hint,
1230 std::make_pair(std::make_pair(cell->level(),
1232 std::move(particle)));
1244 typename std::multimap<internal::LevelInd,
1246 position_hint = particles.end();
1247 for (
auto &particle : loaded_particles_on_cell)
1250 mapping->transform_real_to_unit_cell(cell,
1251 particle.get_location());
1252 particle.set_reference_location(p_unit);
1257 # ifdef DEAL_II_WITH_CXX14 1259 particles.emplace_hint(position_hint,
1260 std::make_pair(cell->level(),
1262 std::move(particle));
1265 particles.insert(position_hint,
1266 std::make_pair(std::make_pair(cell->level(),
1268 std::move(particle)));
1281 typename std::multimap<internal::LevelInd,
1284 for (
unsigned int child_index = 0;
1285 child_index < GeometryInfo<dim>::max_children_per_cell;
1289 child = cell->child(child_index);
1290 position_hints[child_index] = particles.upper_bound(
1291 std::make_pair(child->level(), child->index()));
1294 for (
auto &particle : loaded_particles_on_cell)
1296 for (
unsigned int child_index = 0;
1297 child_index < GeometryInfo<dim>::max_children_per_cell;
1301 child = cell->child(child_index);
1306 mapping->transform_real_to_unit_cell(
1307 child, particle.get_location());
1310 particle.set_reference_location(p_unit);
1316 # ifdef DEAL_II_WITH_CXX14 1317 position_hints[child_index] =
1318 particles.emplace_hint(
1319 position_hints[child_index],
1320 std::make_pair(child->level(), child->index()),
1321 std::move(particle));
1323 position_hints[child_index] = particles.insert(
1324 position_hints[child_index],
1325 std::make_pair(std::make_pair(child->level(),
1327 std::move(particle)));
1332 ++position_hints[child_index];
1352 DEAL_II_NAMESPACE_CLOSE
1354 DEAL_II_NAMESPACE_OPEN
1356 #ifdef DEAL_II_WITH_P4EST 1357 # include "particle_handler.inst" 1360 DEAL_II_NAMESPACE_CLOSE
types::particle_index n_global_particles() const
void sort_particles_into_subdomains_and_cells()
unsigned int register_data_attach(const std::function< std::vector< char >(const cell_iterator &, const CellStatus)> &pack_callback, const bool returns_variable_size_data)
static const unsigned int invalid_unsigned_int
#define AssertDimension(dim1, dim2)
types::particle_index n_locally_owned_particles() const
std::vector< char > store_particles(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const typename Triangulation< dim, spacedim >::CellStatus status) const
void remove_particle(const particle_iterator &particle)
void initialize(const parallel::distributed::Triangulation< dim, spacedim > &tria, const Mapping< dim, spacedim > &mapping, const unsigned int n_properties=0)
types::particle_index n_global_max_particles_per_cell() const
particle_iterator insert_particle(const Particle< dim, spacedim > &particle, const typename Triangulation< dim, spacedim >::active_cell_iterator &cell)
void update_cached_numbers()
const ArrayView< double > get_properties()
#define AssertThrow(cond, exc)
numbers::NumberTraits< Number >::real_type norm() const
std::array< unsigned int, 4 > binary_type
cell_iterator end() const
unsigned int particle_index
void register_load_callback_function(const bool serialization)
static ::ExceptionBase & ExcMessage(std::string arg1)
unsigned int subdomain_id
void insert_particles(const std::multimap< typename Triangulation< dim, spacedim >::active_cell_iterator, Particle< dim, spacedim >> &particles)
particle_iterator begin_ghost() const
unsigned int n_properties_per_particle() const
#define Assert(cond, exc)
particle_iterator_range particles_in_cell(const typename Triangulation< dim, spacedim >::active_cell_iterator &cell)
Abstract base class for mapping classes.
void register_store_callback_function()
void notify_ready_to_unpack(const unsigned int handle, const std::function< void(const cell_iterator &, const CellStatus, const boost::iterator_range< std::vector< char >::const_iterator > &)> &unpack_callback)
particle_iterator begin() const
void exchange_ghost_particles()
size_t pack(const T &object, std::vector< char > &dest_buffer, const bool allow_compression=true)
unsigned int n_mpi_processes(const MPI_Comm &mpi_communicator)
const types::subdomain_id artificial_subdomain_id
types::particle_index get_next_free_particle_index() const
boost::iterator_range< particle_iterator > particle_iterator_range
bool has_properties() const
#define AssertThrowMPI(error_code)
PropertyPool & get_property_pool() const
static ::ExceptionBase & ExcNotImplemented()
void send_recv_particles(const std::map< types::subdomain_id, std::vector< particle_iterator >> &particles_to_send, std::multimap< internal::LevelInd, Particle< dim, spacedim >> &received_particles, const std::map< types::subdomain_id, std::vector< typename Triangulation< dim, spacedim >::active_cell_iterator >> &new_cells_for_particles=std::map< types::subdomain_id, std::vector< typename Triangulation< dim, spacedim >::active_cell_iterator >>())
void load_particles(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const typename Triangulation< dim, spacedim >::CellStatus status, const boost::iterator_range< std::vector< char >::const_iterator > &data_range)
particle_iterator end_ghost() const
particle_iterator end() const
void register_additional_store_load_functions(const std::function< std::size_t()> &size_callback, const std::function< void *(const particle_iterator &, void *)> &store_callback, const std::function< const void *(const particle_iterator &, const void *)> &load_callback)
unsigned int n_particles_in_cell(const typename Triangulation< dim, spacedim >::active_cell_iterator &cell) const
static ::ExceptionBase & ExcInternalError()