Loading [MathJax]/extensions/TeX/AMSsymbols.js
 deal.II version GIT relicensing-2640-gd0584bbf39 2025-02-13 20:00:00+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
particle_handler.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2017 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_particles_particle_handler_h
16#define dealii_particles_particle_handler_h
17
18#include <deal.II/base/config.h>
19
25
27
28#include <deal.II/fe/mapping.h>
29
31
36
37#include <boost/range/iterator_range.hpp>
38#include <boost/serialization/map.hpp>
39#include <boost/signals2.hpp>
40
41
43
44namespace Particles
45{
63 template <int dim, int spacedim = dim>
65 {
66 public:
71
75 using particle_iterator_range = boost::iterator_range<particle_iterator>;
76
82
87
98 const unsigned int n_properties = 0);
99
103 virtual ~ParticleHandler();
104
110 void
113 const unsigned int n_properties = 0);
114
135 void
136 copy_from(const ParticleHandler<dim, spacedim> &particle_handler);
137
141 void
142 clear();
143
149 void
151
166 void
167 reserve(const std::size_t n_particles);
168
177 void
179
184 begin() const;
185
190 begin();
191
196 end() const;
197
202 end();
203
208 begin_ghost() const;
209
214 begin_ghost();
215
220 end_ghost() const;
221
226 end_ghost();
227
234 const;
235
247
259 const;
260
266 void
267 remove_particle(const particle_iterator &particle);
268
274 void
275 remove_particles(const std::vector<particle_iterator> &particles);
276
285 const Particle<dim, spacedim> &particle,
287
305 const Point<spacedim> &position,
306 const Point<dim> &reference_position,
307 const types::particle_index particle_index,
309 const ArrayView<const double> &properties = {});
310
317 void
319 const std::multimap<
322
332 void
333 insert_particles(const std::vector<Point<spacedim>> &positions);
334
395 std::map<unsigned int, IndexSet>
397 const std::vector<Point<spacedim>> &positions,
398 const std::vector<std::vector<BoundingBox<spacedim>>>
399 &global_bounding_boxes,
400 const std::vector<std::vector<double>> &properties = {},
401 const std::vector<types::particle_index> &ids = {});
402
433 std::map<unsigned int, IndexSet>
435 const std::vector<Particle<dim, spacedim>> &particles,
436 const std::vector<std::vector<BoundingBox<spacedim>>>
437 &global_bounding_boxes);
438
482 template <typename VectorType>
483 std::enable_if_t<
484 std::is_convertible_v<VectorType *, Function<spacedim> *> == false>
485 set_particle_positions(const VectorType &input_vector,
486 const bool displace_particles = true);
487
508 void
509 set_particle_positions(const std::vector<Point<spacedim>> &new_positions,
510 const bool displace_particles = true);
511
512
531 void
533 const bool displace_particles = true);
534
562 template <typename VectorType>
563 void
564 get_particle_positions(VectorType &output_vector,
565 const bool add_to_output_vector = false);
566
581 void
582 get_particle_positions(std::vector<Point<spacedim>> &positions,
583 const bool add_to_output_vector = false);
584
611 void
613 const std::function<std::size_t()> &size_callback,
614 const std::function<void *(const particle_iterator &, void *)>
616 const std::function<const void *(const particle_iterator &, const void *)>
618
628 n_global_particles() const;
629
638
645
653
672
676 unsigned int
678
694
699 PropertyPool<dim, spacedim> &
700 get_property_pool() const;
701
716 void
718
724 void
725 exchange_ghost_particles(const bool enable_ghost_cache = false);
726
734 void
736
746 void
748
767 void
769
781 void
783
792 void
793 deserialize();
794
801 template <class Archive>
802 void
803 serialize(Archive &ar, const unsigned int version);
804
817 struct Signals
818 {
832 boost::signals2::signal<void(
833 const typename Particles::ParticleIterator<dim, spacedim> &particle,
835 &cell)>
837 };
838
844
845 private:
854 const void *&data,
856
866
876 void
878
885
892
901 std::unique_ptr<PropertyPool<dim, spacedim>> property_pool;
902
907
913 const typename particle_container::iterator owned_particles_end;
914
920 std::vector<typename particle_container::iterator> cells_to_particle_cache;
921
927
933
943
949
961 std::function<std::size_t()> size_callback;
962
973 std::function<void *(const particle_iterator &, void *)> store_callback;
974
985 std::function<const void *(const particle_iterator &, const void *)>
987
995
999 const double tolerance_inside_cell = 1e-12;
1000
1010 std::unique_ptr<GridTools::Cache<dim, spacedim>> triangulation_cache;
1011
1012#ifdef DEAL_II_WITH_MPI
1036 void
1038 const std::map<types::subdomain_id, std::vector<particle_iterator>>
1039 &particles_to_send,
1040 const std::map<
1042 std::vector<
1044 &new_cells_for_particles = std::map<
1046 std::vector<
1048 const bool enable_cache = false);
1049
1061 void
1063 const std::map<types::subdomain_id, std::vector<particle_iterator>>
1064 &particles_to_send);
1065
1066#endif
1067
1075
1080 void
1082
1087 std::vector<boost::signals2::connection> tria_listeners;
1088
1097 void
1099
1105 void
1107
1115 void
1116 notify_ready_to_unpack(const bool serialization);
1117
1123 std::vector<char>
1126 const CellStatus status) const;
1127
1132 void
1135 const CellStatus status,
1136 const boost::iterator_range<std::vector<char>::const_iterator>
1137 &data_range);
1138
1143 typename particle_container::iterator
1145
1150 typename particle_container::iterator
1152
1157 typename particle_container::iterator
1159
1164 typename particle_container::iterator
1166 };
1167
1168
1169
1170 /* ---------------------- inline and template functions ------------------
1171 */
1172
1173 template <int dim, int spacedim>
1176 {
1177 return (const_cast<ParticleHandler<dim, spacedim> *>(this))->begin();
1178 }
1179
1180
1181
1182 template <int dim, int spacedim>
1185 {
1186 return particle_iterator(particle_container_owned_begin(),
1187 *property_pool,
1188 0);
1189 }
1190
1191
1192
1193 template <int dim, int spacedim>
1196 {
1197 return (const_cast<ParticleHandler<dim, spacedim> *>(this))->end();
1198 }
1199
1200
1201
1202 template <int dim, int spacedim>
1205 {
1206 return particle_iterator(particle_container_owned_end(), *property_pool, 0);
1207 }
1208
1209
1210
1211 template <int dim, int spacedim>
1214 {
1215 return (const_cast<ParticleHandler<dim, spacedim> *>(this))->begin_ghost();
1216 }
1217
1218
1219
1220 template <int dim, int spacedim>
1223 {
1224 return particle_iterator(particle_container_ghost_begin(),
1225 *property_pool,
1226 0);
1227 }
1228
1229
1230
1231 template <int dim, int spacedim>
1234 {
1235 return (const_cast<ParticleHandler<dim, spacedim> *>(this))->end_ghost();
1236 }
1237
1238
1239
1240 template <int dim, int spacedim>
1243 {
1244 return particle_iterator(particle_container_ghost_end(), *property_pool, 0);
1245 }
1246
1247
1248
1249 template <int dim, int spacedim>
1252 {
1253 // We should always have at least the three anchor entries in the list of
1254 // particles
1255 Assert(!particles.empty(), ExcInternalError());
1256 typename particle_container::iterator begin =
1257 const_cast<particle_container &>(particles).begin();
1258 return ++begin;
1259 }
1260
1261
1262
1263 template <int dim, int spacedim>
1266 {
1267 Assert(!particles.empty(), ExcInternalError());
1268 return owned_particles_end;
1269 }
1270
1271
1272
1273 template <int dim, int spacedim>
1276 {
1277 Assert(!particles.empty(), ExcInternalError());
1278 typename particle_container::iterator begin = owned_particles_end;
1279 return ++begin;
1280 }
1281
1282
1283
1284 template <int dim, int spacedim>
1287 {
1288 Assert(!particles.empty(), ExcInternalError());
1289 typename particle_container::iterator end =
1290 const_cast<particle_container &>(particles).end();
1291 return --end;
1292 }
1293
1294
1295
1296 template <int dim, int spacedim>
1297 template <class Archive>
1298 inline void
1299 ParticleHandler<dim, spacedim>::serialize(Archive &ar, const unsigned int)
1300 {
1301 // Note that we do not serialize the particle data itself (i.e., the
1302 // 'particles' member variable). Instead we
1303 // use the serialization functionality of the triangulation class, because
1304 // this guarantees that data is immediately shipped to new processes if
1305 // the domain is distributed differently after resuming from a checkpoint.
1306 //
1307 // See @ref step_83 "step-83" for how to serialize ParticleHandler objects.
1308 ar &global_number_of_particles &global_max_particles_per_cell
1309 &next_free_particle_index;
1310 }
1311
1312
1313
1314 template <int dim, int spacedim>
1315 template <typename VectorType>
1316 inline std::enable_if_t<
1317 std::is_convertible_v<VectorType *, Function<spacedim> *> == false>
1319 const VectorType &input_vector,
1320 const bool displace_particles)
1321 {
1322 AssertDimension(input_vector.size(),
1323 get_next_free_particle_index() * spacedim);
1324 for (auto &p : *this)
1325 {
1326 Point<spacedim> &position = p.get_location();
1327 const auto id = p.get_id();
1328
1329 if (displace_particles)
1330 for (unsigned int i = 0; i < spacedim; ++i)
1331 position[i] += input_vector[id * spacedim + i];
1332 else
1333 for (unsigned int i = 0; i < spacedim; ++i)
1334 position[i] = input_vector[id * spacedim + i];
1335 }
1336 sort_particles_into_subdomains_and_cells();
1337 }
1338
1339
1340
1341 template <int dim, int spacedim>
1342 template <typename VectorType>
1343 inline void
1345 VectorType &output_vector,
1346 const bool add_to_output_vector)
1347 {
1348 AssertDimension(output_vector.size(),
1349 get_next_free_particle_index() * spacedim);
1350 for (const auto &p : *this)
1351 {
1352 auto point = p.get_location();
1353 const auto id = p.get_id();
1354 if (add_to_output_vector)
1355 for (unsigned int i = 0; i < spacedim; ++i)
1356 output_vector[id * spacedim + i] += point[i];
1357 else
1358 for (unsigned int i = 0; i < spacedim; ++i)
1359 output_vector[id * spacedim + i] = point[i];
1360 }
1361 if (add_to_output_vector)
1362 output_vector.compress(VectorOperation::add);
1363 else
1364 output_vector.compress(VectorOperation::insert);
1365 }
1366
1367} // namespace Particles
1368
1370
1371#endif
CellStatus
Definition cell_status.h:31
Abstract base class for mapping classes.
Definition mapping.h:320
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)
types::particle_index n_global_particles() const
particle_iterator end_ghost() const
unsigned int global_max_particles_per_cell
internal::GhostParticlePartitioner< dim, spacedim > ghost_particles_cache
particle_container::iterator particle_container_ghost_begin() const
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
particle_container::iterator particle_container_ghost_end() const
boost::iterator_range< particle_iterator > particle_iterator_range
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)
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
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 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)
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
IndexSet locally_owned_particle_ids() const
Definition point.h:113
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:518
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:519
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcInternalError()
std::vector< index_type > data
Definition mpi.cc:740
boost::signals2::signal< void(const typename Particles::ParticleIterator< dim, spacedim > &particle, const typename Triangulation< dim, spacedim >::active_cell_iterator &cell)> particle_lost