31 template <
int dim,
int spacedim>
33 pack_particles(std::vector<Particle<dim, spacedim>> &particles)
35 std::vector<char> buffer;
37 if (particles.size() == 0)
40 buffer.resize(particles.size() *
41 particles.front().serialized_size_in_bytes());
42 void *current_data = buffer.data();
44 for (
const auto &particle : particles)
46 particle.write_data(current_data);
52 template <
int dim,
int spacedim>
53 std::vector<Particle<dim, spacedim>>
55 const boost::iterator_range<std::vector<char>::const_iterator>
57 PropertyPool &property_pool)
59 std::vector<Particle<dim, spacedim>> particles;
61 if (data_range.empty())
64 Particle<dim, spacedim> particle;
65 particle.set_property_pool(property_pool);
66 const unsigned int particle_size = particle.serialized_size_in_bytes();
68 particles.reserve(data_range.size() / particle_size);
70 const void *data =
static_cast<const void *
>(&(*data_range.begin()));
72 while (data < &(*data_range.end()))
74 particles.emplace_back(data, &property_pool);
78 data == &(*data_range.end()),
80 "The particle data could not be deserialized successfully. "
81 "Check that when deserializing the particles you expect the same "
82 "number of properties that were serialized."));
88 template <
int dim,
int spacedim>
93 , global_number_of_particles(0)
94 , global_max_particles_per_cell(0)
95 , next_free_particle_index(0)
105 template <
int dim,
int spacedim>
109 const unsigned int n_properties)
111 , mapping(&mapping, typeid(*this).name())
114 , global_number_of_particles(0)
115 , global_max_particles_per_cell(0)
116 , next_free_particle_index(0)
124 std_cxx14::make_unique<GridTools::Cache<dim, spacedim>>(
triangulation,
130 template <
int dim,
int spacedim>
135 const unsigned int n_properties)
138 mapping = &new_mapping;
141 property_pool = std_cxx14::make_unique<PropertyPool>(n_properties);
145 triangulation_cache =
146 std_cxx14::make_unique<GridTools::Cache<dim, spacedim>>(new_triangulation,
152 template <
int dim,
int spacedim>
157 global_number_of_particles = 0;
158 next_free_particle_index = 0;
159 global_max_particles_per_cell = 0;
164 template <
int dim,
int spacedim>
173 template <
int dim,
int spacedim>
178 unsigned int local_max_particles_per_cell = 0;
179 unsigned int current_particles_per_cell = 0;
185 locally_highest_index =
186 std::max(locally_highest_index, particle->get_id());
188 if (particle->get_surrounding_cell(*
triangulation) != current_cell)
190 current_particles_per_cell = 0;
191 current_cell = particle->get_surrounding_cell(*
triangulation);
194 ++current_particles_per_cell;
195 local_max_particles_per_cell =
196 std::max(local_max_particles_per_cell, current_particles_per_cell);
199 if (
const auto parallel_triangulation =
204 particles.size(), parallel_triangulation->get_communicator());
205 next_free_particle_index =
206 global_number_of_particles == 0 ?
209 locally_highest_index,
210 parallel_triangulation->get_communicator()) +
213 local_max_particles_per_cell,
214 parallel_triangulation->get_communicator());
218 global_number_of_particles = particles.size();
219 next_free_particle_index =
220 global_number_of_particles == 0 ? 0 : locally_highest_index + 1;
221 global_max_particles_per_cell = local_max_particles_per_cell;
227 template <
int dim,
int spacedim>
236 template <
int dim,
int spacedim>
240 return particle_iterator(particles, particles.begin());
245 template <
int dim,
int spacedim>
254 template <
int dim,
int spacedim>
258 return particle_iterator(particles, particles.end());
263 template <
int dim,
int spacedim>
272 template <
int dim,
int spacedim>
276 return particle_iterator(ghost_particles, ghost_particles.begin());
281 template <
int dim,
int spacedim>
290 template <
int dim,
int spacedim>
294 return particle_iterator(ghost_particles, ghost_particles.end());
299 template <
int dim,
int spacedim>
306 ->particles_in_cell(cell);
311 template <
int dim,
int spacedim>
317 std::make_pair(cell->level(), cell->index());
319 if (cell->is_ghost())
321 const auto particles_in_cell = ghost_particles.equal_range(level_index);
327 const auto particles_in_cell = particles.equal_range(level_index);
335 template <
int dim,
int spacedim>
340 particles.erase(particle->particle);
345 template <
int dim,
int spacedim>
358 particle_it->set_property_pool(*property_pool);
361 for (
unsigned int n = 0; n < particle.
get_properties().size(); ++n)
369 template <
int dim,
int spacedim>
376 for (
auto particle = new_particles.begin(); particle != new_particles.end();
381 particle->first->index()),
384 update_cached_numbers();
389 template <
int dim,
int spacedim>
394 update_cached_numbers();
401 get_next_free_particle_index();
404 #ifdef DEAL_II_WITH_MPI
405 if (
const auto parallel_triangulation =
410 const int ierr = MPI_Scan(&particles_to_add_locally,
415 parallel_triangulation->get_communicator());
417 local_start_index -= particles_to_add_locally;
421 local_start_index += local_next_particle_index;
426 auto &cells = std::get<0>(point_locations);
427 auto &local_positions = std::get<1>(point_locations);
428 auto &index_map = std::get<2>(point_locations);
430 if (cells.size() == 0)
434 particles.find(std::make_pair(cells[0]->
level(), cells[0]->index()));
435 for (
unsigned int i = 0; i < cells.size(); ++i)
438 for (
unsigned int p = 0; p < local_positions[i].size(); ++p)
440 hint = particles.insert(
442 std::make_pair(current_cell,
444 local_positions[i][p],
450 update_cached_numbers();
455 template <
int dim,
int spacedim>
456 std::map<unsigned int, IndexSet>
460 & global_bounding_boxes,
461 const std::vector<std::vector<double>> &properties)
463 if (!properties.empty())
467 for (
const auto &p : properties)
476 (tria !=
nullptr ? tria->get_communicator() : MPI_COMM_WORLD);
483 const auto n_global_properties =
487 const auto n_particles_per_proc =
493 unsigned int particle_start_index = 0;
494 for (
unsigned int process = 0; process < particle_start_indices.size();
497 particle_start_indices[process] = particle_start_index;
498 particle_start_index += n_particles_per_proc[process];
502 const auto cells_positions_and_index_maps =
505 global_bounding_boxes);
509 const auto &local_cells_containing_particles =
510 std::get<0>(cells_positions_and_index_maps);
514 const auto &local_reference_positions =
515 std::get<1>(cells_positions_and_index_maps);
518 const auto &original_indices_of_local_particles =
519 std::get<2>(cells_positions_and_index_maps);
522 const auto &local_positions = std::get<3>(cells_positions_and_index_maps);
524 const auto &calling_process_indices =
525 std::get<4>(cells_positions_and_index_maps);
529 std::map<unsigned int, IndexSet> original_process_to_local_particle_indices;
530 for (
unsigned int i_cell = 0;
531 i_cell < local_cells_containing_particles.size();
534 for (
unsigned int i_particle = 0;
535 i_particle < local_positions[i_cell].size();
538 const unsigned int local_id_on_calling_process =
539 original_indices_of_local_particles[i_cell][i_particle];
540 const unsigned int calling_process =
541 calling_process_indices[i_cell][i_particle];
543 if (original_process_to_local_particle_indices.find(
545 original_process_to_local_particle_indices.end())
546 original_process_to_local_particle_indices.insert(
548 IndexSet(n_particles_per_proc[calling_process])});
550 original_process_to_local_particle_indices[calling_process]
551 .add_index(local_id_on_calling_process);
554 for (
auto &process_and_particle_indices :
555 original_process_to_local_particle_indices)
556 process_and_particle_indices.second.compress();
564 std::map<unsigned int, std::vector<std::vector<double>>>
565 locally_owned_properties_from_other_processes;
567 if (n_global_properties > 0)
572 comm, original_process_to_local_particle_indices);
574 std::map<unsigned int, std::vector<std::vector<double>>>
575 non_locally_owned_properties;
578 for (
const auto &it : send_to_cpu)
580 std::vector<std::vector<double>> properties_to_send(
581 it.second.n_elements(),
582 std::vector<double>(n_properties_per_particle()));
583 unsigned int index = 0;
584 for (
const auto &el : it.second)
585 properties_to_send[index++] = properties[el];
586 non_locally_owned_properties.insert({it.first, properties_to_send});
591 locally_owned_properties_from_other_processes =
595 original_process_to_local_particle_indices.size());
600 std::multimap<typename Triangulation<dim, spacedim>::active_cell_iterator,
605 for (
unsigned int i_cell = 0;
606 i_cell < local_cells_containing_particles.size();
609 for (
unsigned int i_particle = 0;
610 i_particle < local_positions[i_cell].size();
613 const unsigned int local_id_on_calling_process =
614 original_indices_of_local_particles[i_cell][i_particle];
615 const unsigned int calling_process =
616 calling_process_indices[i_cell][i_particle];
618 const unsigned int particle_id =
619 local_id_on_calling_process +
620 particle_start_indices[calling_process];
623 local_positions[i_cell][i_particle],
624 local_reference_positions[i_cell][i_particle],
627 if (n_global_properties > 0)
629 const unsigned int index_within_set =
630 original_process_to_local_particle_indices[calling_process]
631 .index_within_set(local_id_on_calling_process);
633 const auto &this_particle_properties =
634 locally_owned_properties_from_other_processes
635 [calling_process][index_within_set];
641 particles.emplace(local_cells_containing_particles[i_cell],
646 this->insert_particles(particles);
648 return original_process_to_local_particle_indices;
653 template <
int dim,
int spacedim>
657 return global_number_of_particles;
662 template <
int dim,
int spacedim>
666 return global_max_particles_per_cell;
671 template <
int dim,
int spacedim>
675 return particles.size();
680 template <
int dim,
int spacedim>
684 return property_pool->n_properties_per_slot();
689 template <
int dim,
int spacedim>
696 std::make_pair(cell->level(), cell->index());
698 if (cell->is_locally_owned())
699 return particles.count(found_cell);
700 else if (cell->is_ghost())
701 return ghost_particles.count(found_cell);
702 else if (cell->is_artificial())
710 template <
int dim,
int spacedim>
714 return next_free_particle_index;
719 template <
int dim,
int spacedim>
724 for (
const auto &p : *
this)
725 set.add_index(p.get_id());
732 template <
int dim,
int spacedim>
736 const bool add_to_output_vector)
742 for (
auto it =
begin(); it !=
end(); ++it, ++i)
744 if (add_to_output_vector)
745 positions[i] = positions[i] + it->get_location();
747 positions[i] = it->get_location();
753 template <
int dim,
int spacedim>
757 const bool displace_particles)
763 for (
auto it =
begin(); it !=
end(); ++it, ++i)
765 if (displace_particles)
766 it->set_location(it->get_location() + new_positions[i]);
768 it->set_location(new_positions[i]);
770 sort_particles_into_subdomains_and_cells();
775 template <
int dim,
int spacedim>
779 const bool displace_particles)
786 for (
auto &particle : *
this)
789 function.vector_value(particle_location, new_position);
790 if (displace_particles)
791 for (
unsigned int d = 0;
d < spacedim; ++
d)
792 particle_location[
d] += new_position[
d];
794 for (
unsigned int d = 0;
d < spacedim; ++
d)
795 particle_location[
d] = new_position[
d];
796 particle.set_location(particle_location);
798 sort_particles_into_subdomains_and_cells();
803 template <
int dim,
int spacedim>
807 return *property_pool;
825 compare_particle_association(
826 const unsigned int a,
827 const unsigned int b,
831 const double scalar_product_a = center_directions[a] * particle_direction;
832 const double scalar_product_b = center_directions[
b] * particle_direction;
837 return (scalar_product_a > scalar_product_b);
843 template <
int dim,
int spacedim>
855 std::vector<particle_iterator> particles_out_of_cell;
856 particles_out_of_cell.reserve(n_locally_owned_particles());
867 mapping->transform_real_to_unit_cell(cell, it->get_location());
870 it->set_reference_location(p_unit);
875 particles_out_of_cell.push_back(it);
881 particles_out_of_cell.push_back(it);
891 std::vector<std::pair<internal::LevelInd, Particle<dim, spacedim>>>
893 std::map<types::subdomain_id, std::vector<particle_iterator>>
897 std::vector<typename Triangulation<dim, spacedim>::active_cell_iterator>>
906 sorted_particles.reserve(
907 static_cast<vector_size
>(particles_out_of_cell.size() * 1.25));
909 std::set<types::subdomain_id> ghost_owners;
910 if (
const auto parallel_triangulation =
913 ghost_owners = parallel_triangulation->ghost_owners();
915 for (
const auto ghost_owner : ghost_owners)
916 moved_particles[ghost_owner].reserve(
917 static_cast<vector_size
>(particles_out_of_cell.size() * 0.25));
918 for (
const auto ghost_owner : ghost_owners)
919 moved_cells[ghost_owner].reserve(
920 static_cast<vector_size
>(particles_out_of_cell.size() * 0.25));
925 std::set<typename Triangulation<dim, spacedim>::active_cell_iterator>>
926 vertex_to_cells = triangulation_cache->get_vertex_to_cell_map();
930 std::vector<std::vector<Tensor<1, spacedim>>> vertex_to_cell_centers =
931 triangulation_cache->get_vertex_to_cell_centers_directions();
933 std::vector<unsigned int> neighbor_permutation;
936 typename std::vector<particle_iterator>::iterator
937 it = particles_out_of_cell.begin(),
938 end_particle = particles_out_of_cell.end();
940 for (; it != end_particle; ++it)
944 bool found_cell =
false;
951 const unsigned int closest_vertex =
952 GridTools::find_closest_vertex_of_cell<dim, spacedim>(
953 current_cell, (*it)->get_location());
955 (*it)->get_location() - current_cell->vertex(closest_vertex);
956 vertex_to_particle /= vertex_to_particle.
norm();
958 const unsigned int closest_vertex_index =
959 current_cell->vertex_index(closest_vertex);
960 const unsigned int n_neighbor_cells =
961 vertex_to_cells[closest_vertex_index].size();
963 neighbor_permutation.resize(n_neighbor_cells);
964 for (
unsigned int i = 0; i < n_neighbor_cells; ++i)
965 neighbor_permutation[i] = i;
967 const auto cell_centers =
968 vertex_to_cell_centers[closest_vertex_index];
969 std::sort(neighbor_permutation.begin(),
970 neighbor_permutation.end(),
971 [&vertex_to_particle, &cell_centers](
const unsigned int a,
972 const unsigned int b) {
973 return compare_particle_association(a,
981 for (
unsigned int i = 0; i < n_neighbor_cells; ++i)
985 typename std::set<typename Triangulation<dim, spacedim>::
986 active_cell_iterator>::const_iterator
987 cell = vertex_to_cells[closest_vertex_index].begin();
991 mapping->transform_real_to_unit_cell(*cell,
992 (*it)->get_location());
995 current_cell = *cell;
996 current_reference_position = p_unit;
1012 const std::pair<const typename Triangulation<dim, spacedim>::
1013 active_cell_iterator,
1015 current_cell_and_position =
1016 GridTools::find_active_cell_around_point<>(
1018 current_cell = current_cell_and_position.first;
1019 current_reference_position = current_cell_and_position.second;
1021 catch (GridTools::ExcPointNotFound<spacedim> &)
1031 (*it)->set_reference_location(current_reference_position);
1035 if (current_cell->is_locally_owned())
1037 sorted_particles.push_back(
1039 current_cell->index()),
1040 (*it)->particle->second));
1044 moved_particles[current_cell->subdomain_id()].push_back(*it);
1045 moved_cells[current_cell->subdomain_id()].push_back(current_cell);
1052 std::multimap<internal::LevelInd, Particle<dim, spacedim>>
1053 sorted_particles_map;
1056 #ifdef DEAL_II_WITH_MPI
1057 if (
const auto parallel_triangulation =
1062 parallel_triangulation->get_communicator()) > 1)
1063 send_recv_particles(moved_particles,
1064 sorted_particles_map,
1069 sorted_particles_map.insert(sorted_particles.begin(),
1070 sorted_particles.end());
1072 for (
unsigned int i = 0; i < particles_out_of_cell.size(); ++i)
1073 remove_particle(particles_out_of_cell[i]);
1075 particles.insert(sorted_particles_map.begin(), sorted_particles_map.end());
1076 update_cached_numbers();
1081 template <
int dim,
int spacedim>
1086 const auto parallel_triangulation =
1089 if (parallel_triangulation !=
nullptr)
1092 parallel_triangulation->get_communicator()) == 1)
1098 #ifdef DEAL_II_WITH_MPI
1100 ghost_particles.
clear();
1102 std::map<types::subdomain_id, std::vector<particle_iterator>>
1103 ghost_particles_by_domain;
1105 const std::set<types::subdomain_id> ghost_owners =
1106 parallel_triangulation->ghost_owners();
1107 for (
const auto ghost_owner : ghost_owners)
1108 ghost_particles_by_domain[ghost_owner].reserve(
1110 particles.size() * 0.25));
1112 std::vector<std::set<unsigned int>> vertex_to_neighbor_subdomain(
1115 for (
const auto &cell :
triangulation->active_cell_iterators())
1117 if (cell->is_ghost())
1119 vertex_to_neighbor_subdomain[cell->vertex_index(v)].insert(
1120 cell->subdomain_id());
1123 for (
const auto &cell :
triangulation->active_cell_iterators())
1125 if (!cell->is_ghost())
1127 std::set<unsigned int> cell_to_neighbor_subdomain;
1130 cell_to_neighbor_subdomain.insert(
1131 vertex_to_neighbor_subdomain[cell->vertex_index(v)].begin(),
1132 vertex_to_neighbor_subdomain[cell->vertex_index(v)].end());
1135 if (cell_to_neighbor_subdomain.size() > 0)
1138 particles_in_cell(cell);
1140 for (
const auto domain : cell_to_neighbor_subdomain)
1142 for (
typename particle_iterator_range::iterator particle =
1143 particle_range.begin();
1144 particle != particle_range.end();
1146 ghost_particles_by_domain[domain].push_back(particle);
1152 send_recv_particles(ghost_particles_by_domain, ghost_particles);
1158 #ifdef DEAL_II_WITH_MPI
1159 template <
int dim,
int spacedim>
1165 &received_particles,
1171 const auto parallel_triangulation =
1175 parallel_triangulation,
1177 "This function is only implemented for parallel::Triangulation objects."));
1180 const std::set<types::subdomain_id> ghost_owners =
1181 parallel_triangulation->ghost_owners();
1182 const std::vector<types::subdomain_id> neighbors(ghost_owners.begin(),
1183 ghost_owners.end());
1184 const unsigned int n_neighbors = neighbors.size();
1186 if (send_cells.size() != 0)
1192 particles_to_send.end(),
1197 for (
auto send_particles = particles_to_send.begin();
1198 send_particles != particles_to_send.end();
1200 Assert(ghost_owners.find(send_particles->first) != ghost_owners.end(),
1203 unsigned int n_send_particles = 0;
1204 for (
auto send_particles = particles_to_send.begin();
1205 send_particles != particles_to_send.end();
1207 n_send_particles += send_particles->second.size();
1213 std::vector<unsigned int> n_send_data(n_neighbors, 0);
1214 std::vector<unsigned int> send_offsets(n_neighbors, 0);
1215 std::vector<char> send_data;
1220 if (n_send_particles > 0)
1223 const unsigned int particle_size =
1224 begin()->serialized_size_in_bytes() + cellid_size +
1225 (size_callback ? size_callback() : 0);
1226 send_data.resize(n_send_particles * particle_size);
1227 void *data =
static_cast<void *
>(&send_data.front());
1230 for (
unsigned int i = 0; i < n_neighbors; ++i)
1232 send_offsets[i] =
reinterpret_cast<std::size_t
>(data) -
1233 reinterpret_cast<std::size_t
>(&send_data.front());
1235 for (
unsigned int j = 0;
1236 j < particles_to_send.at(neighbors[i]).size();
1242 if (send_cells.size() == 0)
1244 particles_to_send.at(neighbors[i])[j]->get_surrounding_cell(
1247 cell = send_cells.at(neighbors[i])[j];
1250 cell->id().template to_binary<dim>();
1251 memcpy(data, &cellid, cellid_size);
1252 data =
static_cast<char *
>(data) + cellid_size;
1254 particles_to_send.at(neighbors[i])[j]->write_data(data);
1257 store_callback(particles_to_send.at(neighbors[i])[j], data);
1259 n_send_data[i] =
reinterpret_cast<std::size_t
>(data) -
1261 reinterpret_cast<std::size_t
>(&send_data.front());
1266 std::vector<unsigned int> n_recv_data(n_neighbors);
1267 std::vector<unsigned int> recv_offsets(n_neighbors);
1274 std::vector<MPI_Request> n_requests(2 * n_neighbors);
1275 for (
unsigned int i = 0; i < n_neighbors; ++i)
1277 const int ierr = MPI_Irecv(&(n_recv_data[i]),
1282 parallel_triangulation->get_communicator(),
1283 &(n_requests[2 * i]));
1286 for (
unsigned int i = 0; i < n_neighbors; ++i)
1288 const int ierr = MPI_Isend(&(n_send_data[i]),
1293 parallel_triangulation->get_communicator(),
1294 &(n_requests[2 * i + 1]));
1298 MPI_Waitall(2 * n_neighbors, n_requests.data(), MPI_STATUSES_IGNORE);
1303 unsigned int total_recv_data = 0;
1304 for (
unsigned int neighbor_id = 0; neighbor_id < n_neighbors; ++neighbor_id)
1306 recv_offsets[neighbor_id] = total_recv_data;
1307 total_recv_data += n_recv_data[neighbor_id];
1311 std::vector<char> recv_data(total_recv_data);
1315 std::vector<MPI_Request> requests(2 * n_neighbors);
1316 unsigned int send_ops = 0;
1317 unsigned int recv_ops = 0;
1322 for (
unsigned int i = 0; i < n_neighbors; ++i)
1323 if (n_recv_data[i] > 0)
1326 MPI_Irecv(&(recv_data[recv_offsets[i]]),
1331 parallel_triangulation->get_communicator(),
1332 &(requests[send_ops]));
1337 for (
unsigned int i = 0; i < n_neighbors; ++i)
1338 if (n_send_data[i] > 0)
1341 MPI_Isend(&(send_data[send_offsets[i]]),
1346 parallel_triangulation->get_communicator(),
1347 &(requests[send_ops + recv_ops]));
1352 MPI_Waitall(send_ops + recv_ops, requests.data(), MPI_STATUSES_IGNORE);
1358 const void *recv_data_it =
static_cast<const void *
>(recv_data.data());
1360 while (
reinterpret_cast<std::size_t
>(recv_data_it) -
1361 reinterpret_cast<std::size_t
>(recv_data.data()) <
1365 memcpy(&binary_cellid, recv_data_it, cellid_size);
1366 const CellId id(binary_cellid);
1367 recv_data_it =
static_cast<const char *
>(recv_data_it) + cellid_size;
1374 recv_particle = received_particles.insert(std::make_pair(
1384 AssertThrow(recv_data_it == recv_data.data() + recv_data.size(),
1386 "The amount of data that was read into new particles "
1387 "does not match the amount of data sent around."));
1393 template <
int dim,
int spacedim>
1396 const std::function<std::size_t()> & size_callb,
1401 size_callback = size_callb;
1402 store_callback = store_callb;
1403 load_callback = load_callb;
1408 template <
int dim,
int spacedim>
1413 *non_const_triangulation =
1416 *
>(&(*triangulation)));
1417 (void)non_const_triangulation;
1421 #ifdef DEAL_II_WITH_P4EST
1425 update_cached_numbers();
1427 if (global_max_particles_per_cell > 0)
1429 const auto callback_function =
1434 return this->store_particles(cell_iterator, cell_status);
1438 callback_function,
true);
1445 template <
int dim,
int spacedim>
1448 const bool serialization)
1455 *non_const_triangulation =
1458 *
>(&(*triangulation)));
1459 (void)non_const_triangulation;
1463 #ifdef DEAL_II_WITH_P4EST
1468 if (serialization && (global_max_particles_per_cell > 0))
1470 const auto callback_function =
1475 return this->store_particles(cell_iterator, cell_status);
1479 callback_function,
true);
1485 const auto callback_function =
1490 const boost::iterator_range<std::vector<char>::const_iterator>
1492 this->load_particles(cell_iterator, cell_status, range_iterator);
1501 update_cached_numbers();
1504 (void)serialization;
1510 template <
int dim,
int spacedim>
1516 std::vector<Particle<dim, spacedim>> stored_particles_on_cell;
1525 unsigned int n_particles = 0;
1529 const auto particles_in_cell =
1530 (cell->is_ghost() ? ghost_particles.equal_range(level_index) :
1531 particles.equal_range(level_index));
1533 n_particles = n_particles_in_cell(cell);
1534 stored_particles_on_cell.reserve(n_particles);
1537 particles_in_cell.first,
1538 particles_in_cell.second,
1539 [&stored_particles_on_cell](
1542 stored_particles_on_cell.push_back(particle.second);
1553 unsigned int n_particles = 0;
1555 for (
unsigned int child_index = 0;
1556 child_index < GeometryInfo<dim>::max_children_per_cell;
1560 child = cell->child(child_index);
1561 n_particles += n_particles_in_cell(child);
1564 stored_particles_on_cell.reserve(n_particles);
1566 for (
unsigned int child_index = 0;
1567 child_index < GeometryInfo<dim>::max_children_per_cell;
1571 child = cell->child(child_index);
1574 const auto particles_in_cell =
1575 (child->is_ghost() ?
1576 ghost_particles.equal_range(level_index) :
1577 particles.equal_range(level_index));
1580 particles_in_cell.first,
1581 particles_in_cell.second,
1582 [&stored_particles_on_cell](
1585 stored_particles_on_cell.push_back(particle.second);
1598 return pack_particles(stored_particles_on_cell);
1601 template <
int dim,
int spacedim>
1606 const boost::iterator_range<std::vector<char>::const_iterator> &data_range)
1610 std::vector<Particle<dim, spacedim>> loaded_particles_on_cell =
1611 unpack_particles<dim, spacedim>(data_range, *property_pool);
1616 for (
auto &particle : loaded_particles_on_cell)
1617 particle.set_property_pool(*property_pool);
1623 auto position_hint = particles.end();
1624 for (
const auto &particle : loaded_particles_on_cell)
1630 #ifdef DEAL_II_WITH_CXX14
1632 particles.emplace_hint(position_hint,
1633 std::make_pair(cell->level(),
1635 std::move(particle));
1638 particles.insert(position_hint,
1639 std::make_pair(std::make_pair(cell->level(),
1641 std::move(particle)));
1655 position_hint = particles.end();
1656 for (
auto &particle : loaded_particles_on_cell)
1659 mapping->transform_real_to_unit_cell(cell,
1660 particle.get_location());
1661 particle.set_reference_location(p_unit);
1666 #ifdef DEAL_II_WITH_CXX14
1668 particles.emplace_hint(position_hint,
1669 std::make_pair(cell->level(),
1671 std::move(particle));
1674 particles.insert(position_hint,
1675 std::make_pair(std::make_pair(cell->level(),
1677 std::move(particle)));
1693 for (
unsigned int child_index = 0;
1694 child_index < GeometryInfo<dim>::max_children_per_cell;
1698 child = cell->child(child_index);
1699 position_hints[child_index] = particles.upper_bound(
1700 std::make_pair(child->level(), child->index()));
1703 for (
auto &particle : loaded_particles_on_cell)
1705 for (
unsigned int child_index = 0;
1706 child_index < GeometryInfo<dim>::max_children_per_cell;
1710 child = cell->child(child_index);
1715 mapping->transform_real_to_unit_cell(
1716 child, particle.get_location());
1719 particle.set_reference_location(p_unit);
1725 #ifdef DEAL_II_WITH_CXX14
1726 position_hints[child_index] =
1727 particles.emplace_hint(
1728 position_hints[child_index],
1729 std::make_pair(child->level(), child->index()),
1730 std::move(particle));
1732 position_hints[child_index] = particles.insert(
1733 position_hints[child_index],
1734 std::make_pair(std::make_pair(child->level(),
1736 std::move(particle)));
1741 ++position_hints[child_index];
1759 #include "particle_handler.inst"