15#ifndef dealii_particles_particle_handler_h
16#define dealii_particles_particle_handler_h
37#include <boost/range/iterator_range.hpp>
38#include <boost/serialization/map.hpp>
61 template <
int dim,
int spacedim = dim>
96 const unsigned int n_properties = 0);
111 const unsigned int n_properties = 0);
393 std::map<unsigned int, IndexSet>
398 const std::vector<std::vector<double>> &properties = {},
399 const std::vector<types::particle_index> &ids = {});
431 std::map<unsigned int, IndexSet>
480 template <
typename VectorType>
482 std::is_convertible_v<VectorType *, Function<spacedim> *> ==
false>
560 template <
typename VectorType>
799 template <
class Archive>
830 boost::signals2::signal<
void(
1010#ifdef DEAL_II_WITH_MPI
1134 const boost::iterator_range<std::vector<char>::const_iterator>
1141 typename particle_container::iterator
1148 typename particle_container::iterator
1155 typename particle_container::iterator
1162 typename particle_container::iterator
1171 template <
int dim,
int spacedim>
1180 template <
int dim,
int spacedim>
1191 template <
int dim,
int spacedim>
1200 template <
int dim,
int spacedim>
1209 template <
int dim,
int spacedim>
1218 template <
int dim,
int spacedim>
1229 template <
int dim,
int spacedim>
1238 template <
int dim,
int spacedim>
1247 template <
int dim,
int spacedim>
1254 typename particle_container::iterator begin =
1261 template <
int dim,
int spacedim>
1266 return owned_particles_end;
1271 template <
int dim,
int spacedim>
1276 typename particle_container::iterator begin = owned_particles_end;
1282 template <
int dim,
int spacedim>
1287 typename particle_container::iterator end =
1294 template <
int dim,
int spacedim>
1295 template <
class Archive>
1306 ar &global_number_of_particles &global_max_particles_per_cell
1307 &next_free_particle_index;
1312 template <
int dim,
int spacedim>
1313 template <
typename VectorType>
1314 inline std::enable_if_t<
1315 std::is_convertible_v<VectorType *, Function<spacedim> *> ==
false>
1321 get_next_free_particle_index() * spacedim);
1322 for (
auto &p : *
this)
1325 const auto id = p.get_id();
1328 for (
unsigned int i = 0; i < spacedim; ++i)
1331 for (
unsigned int i = 0; i < spacedim; ++i)
1334 sort_particles_into_subdomains_and_cells();
1339 template <
int dim,
int spacedim>
1340 template <
typename VectorType>
1347 get_next_free_particle_index() * spacedim);
1348 for (
const auto &p : *
this)
1350 auto point = p.get_location();
1351 const auto id = p.get_id();
1353 for (
unsigned int i = 0; i < spacedim; ++i)
1356 for (
unsigned int i = 0; i < spacedim; ++i)
std::list< ParticlesInCell > particle_container
std::function< void *(const particle_iterator &, void *)> store_callback
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)
void exchange_ghost_particles(const bool enable_ghost_cache=false)
unsigned int tria_attached_data_index
types::particle_index n_global_particles() const
particle_iterator end_ghost() const
unsigned int global_max_particles_per_cell
void unpack_after_coarsening_and_refinement()
void register_data_attach()
internal::GhostParticlePartitioner< dim, spacedim > ghost_particles_cache
particle_container::iterator particle_container_ghost_begin() const
void prepare_for_serialization()
unsigned int n_properties_per_particle() const
types::particle_index global_number_of_particles
void get_particle_positions(VectorType &output_vector, const bool add_to_output_vector=false)
particle_container::iterator particle_container_owned_end() const
const double tolerance_inside_cell
void update_cached_numbers()
void prepare_for_coarsening_and_refinement()
particle_container::iterator particle_container_ghost_end() const
boost::iterator_range< particle_iterator > particle_iterator_range
void post_mesh_change_action()
void send_recv_particles(const std::map< types::subdomain_id, std::vector< particle_iterator > > &particles_to_send, 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 > >(), const bool enable_cache=false)
ParticleIterator< dim, spacedim > particle_iterator
ObserverPointer< const Mapping< dim, spacedim >, ParticleHandler< dim, spacedim > > mapping
void send_recv_particles_properties_and_location(const std::map< types::subdomain_id, std::vector< particle_iterator > > &particles_to_send)
types::particle_index number_of_locally_owned_particles
particle_iterator begin_ghost() const
PropertyPool< dim, spacedim > & get_property_pool() const
void initialize(const Triangulation< dim, spacedim > &tria, const Mapping< dim, spacedim > &mapping, const unsigned int n_properties=0)
void serialize(Archive &ar, const unsigned int version)
std::unique_ptr< PropertyPool< dim, spacedim > > property_pool
types::particle_index get_max_local_particle_index() const
void reserve(const std::size_t n_particles)
particle_iterator begin() const
std::map< unsigned int, IndexSet > insert_global_particles(const std::vector< Point< spacedim > > &positions, const std::vector< std::vector< BoundingBox< spacedim > > > &global_bounding_boxes, const std::vector< std::vector< double > > &properties={}, const std::vector< types::particle_index > &ids={})
void notify_ready_to_unpack(const bool serialization)
void connect_to_triangulation_signals()
std::enable_if_t< std::is_convertible_v< VectorType *, Function< spacedim > * >==false > set_particle_positions(const VectorType &input_vector, const bool displace_particles=true)
std::function< std::size_t()> size_callback
std::vector< char > pack_callback(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellStatus status) const
void insert_particles(const std::multimap< typename Triangulation< dim, spacedim >::active_cell_iterator, Particle< dim, spacedim > > &particles)
types::particle_index next_free_particle_index
const particle_container::iterator owned_particles_end
particle_iterator end() const
void remove_particle(const particle_iterator &particle)
types::particle_index n_locally_owned_particles() const
void reset_particle_container(particle_container &particles)
types::particle_index n_particles_in_cell(const typename Triangulation< dim, spacedim >::active_cell_iterator &cell) const
void update_ghost_particles()
typename ParticleAccessor< dim, spacedim >::particle_container particle_container
ObserverPointer< const Triangulation< dim, spacedim >, ParticleHandler< dim, spacedim > > triangulation
particle_iterator_range particles_in_cell(const typename Triangulation< dim, spacedim >::active_cell_iterator &cell)
particle_iterator insert_particle(const Particle< dim, spacedim > &particle, const typename Triangulation< dim, spacedim >::active_cell_iterator &cell)
void sort_particles_into_subdomains_and_cells()
void remove_particles(const std::vector< particle_iterator > &particles)
types::particle_index n_global_max_particles_per_cell() const
std::vector< boost::signals2::connection > tria_listeners
particle_container::iterator particle_container_owned_begin() const
void unpack_callback(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellStatus status, const boost::iterator_range< std::vector< char >::const_iterator > &data_range)
virtual ~ParticleHandler()
std::vector< typename particle_container::iterator > cells_to_particle_cache
void copy_from(const ParticleHandler< dim, spacedim > &particle_handler)
std::function< const void *(const particle_iterator &, const void *)> load_callback
std::unique_ptr< GridTools::Cache< dim, spacedim > > triangulation_cache
types::particle_index get_next_free_particle_index() const
particle_container particles
IndexSet locally_owned_particle_ids() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcInternalError()
std::vector< index_type > data
boost::signals2::signal< void(const typename Particles::ParticleIterator< dim, spacedim > &particle, const typename Triangulation< dim, spacedim >::active_cell_iterator &cell)> particle_lost