Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
particle_handler.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2017 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16#ifndef dealii_particles_particle_handler_h
17#define dealii_particles_particle_handler_h
18
19#include <deal.II/base/config.h>
20
26
28
29#include <deal.II/fe/mapping.h>
30
32
37
38#include <boost/range/iterator_range.hpp>
39#include <boost/serialization/map.hpp>
40
42
43namespace Particles
44{
64 template <int dim, int spacedim = dim>
66 {
67 public:
72
76 using particle_iterator_range = boost::iterator_range<particle_iterator>;
77
83
88
99 const unsigned int n_properties = 0);
100
104 virtual ~ParticleHandler();
105
111 void
114 const unsigned int n_properties = 0);
115
136 void
137 copy_from(const ParticleHandler<dim, spacedim> &particle_handler);
138
142 void
143 clear();
144
150 void
152
167 void
168 reserve(const std::size_t n_particles);
169
178 void
180
185 begin() const;
186
191 begin();
192
197 end() const;
198
203 end();
204
209 begin_ghost() const;
210
215 begin_ghost();
216
221 end_ghost() const;
222
227 end_ghost();
228
235 const;
236
248
260 const;
261
267 void
268 remove_particle(const particle_iterator &particle);
269
275 void
276 remove_particles(const std::vector<particle_iterator> &particles);
277
286 const Particle<dim, spacedim> &particle,
288
306 const Point<spacedim> & position,
307 const Point<dim> & reference_position,
308 const types::particle_index particle_index,
310 const ArrayView<const double> &properties = {});
311
318 void
320 const std::multimap<
323
333 void
334 insert_particles(const std::vector<Point<spacedim>> &positions);
335
396 std::map<unsigned int, IndexSet>
398 const std::vector<Point<spacedim>> &positions,
399 const std::vector<std::vector<BoundingBox<spacedim>>>
400 & global_bounding_boxes,
401 const std::vector<std::vector<double>> & properties = {},
402 const std::vector<types::particle_index> &ids = {});
403
434 std::map<unsigned int, IndexSet>
436 const std::vector<Particle<dim, spacedim>> &particles,
437 const std::vector<std::vector<BoundingBox<spacedim>>>
438 &global_bounding_boxes);
439
483 template <class VectorType>
484 std::enable_if_t<
485 std::is_convertible<VectorType *, Function<spacedim> *>::value == false>
486 set_particle_positions(const VectorType &input_vector,
487 const bool displace_particles = true);
488
509 void
510 set_particle_positions(const std::vector<Point<spacedim>> &new_positions,
511 const bool displace_particles = true);
512
513
532 void
534 const bool displace_particles = true);
535
563 template <class VectorType>
564 void
565 get_particle_positions(VectorType &output_vector,
566 const bool add_to_output_vector = false);
567
582 void
583 get_particle_positions(std::vector<Point<spacedim>> &positions,
584 const bool add_to_output_vector = false);
585
612 void
614 const std::function<std::size_t()> &size_callback,
615 const std::function<void *(const particle_iterator &, void *)>
617 const std::function<const void *(const particle_iterator &, const void *)>
619
629 n_global_particles() const;
630
639
646
654
673
677 unsigned int
679
695
700 PropertyPool<dim, spacedim> &
701 get_property_pool() const;
702
717 void
719
725 void
726 exchange_ghost_particles(const bool enable_ghost_cache = false);
727
735 void
737
747 void
749
768 void
770
782 void
784
793 void
794 deserialize();
795
807 void
809
821 void
822 register_load_callback_function(const bool serialization);
823
828 template <class Archive>
829 void
830 serialize(Archive &ar, const unsigned int version);
831
844 struct Signals
845 {
859 boost::signals2::signal<void(
860 const typename Particles::ParticleIterator<dim, spacedim> &particle,
862 &cell)>
864 };
865
871
872 private:
881 const void *& data,
883
892
902 void
904
911
917
926 std::unique_ptr<PropertyPool<dim, spacedim>> property_pool;
927
932
938 const typename particle_container::iterator owned_particles_end;
939
945 std::vector<typename particle_container::iterator> cells_to_particle_cache;
946
952
958
968
974
986 std::function<std::size_t()> size_callback;
987
998 std::function<void *(const particle_iterator &, void *)> store_callback;
999
1010 std::function<const void *(const particle_iterator &, const void *)>
1012
1019 unsigned int handle;
1020
1030 std::unique_ptr<GridTools::Cache<dim, spacedim>> triangulation_cache;
1031
1032#ifdef DEAL_II_WITH_MPI
1056 void
1058 const std::map<types::subdomain_id, std::vector<particle_iterator>>
1059 &particles_to_send,
1060 const std::map<
1062 std::vector<
1064 &new_cells_for_particles = std::map<
1066 std::vector<
1068 const bool enable_cache = false);
1069
1081 void
1083 const std::map<types::subdomain_id, std::vector<particle_iterator>>
1084 &particles_to_send);
1085
1086#endif
1087
1095
1100 void
1102
1107 std::vector<boost::signals2::connection> tria_listeners;
1108
1117 void
1119
1125 void
1127
1135 void
1136 notify_ready_to_unpack(const bool serialization);
1137
1143 std::vector<char>
1146 const typename Triangulation<dim, spacedim>::CellStatus status) const;
1147
1152 void
1155 const typename Triangulation<dim, spacedim>::CellStatus status,
1156 const boost::iterator_range<std::vector<char>::const_iterator>
1157 &data_range);
1158
1163 typename particle_container::iterator
1165
1170 typename particle_container::iterator
1172
1177 typename particle_container::iterator
1179
1184 typename particle_container::iterator
1186 };
1187
1188
1189
1190 /* ---------------------- inline and template functions ------------------
1191 */
1192
1193 template <int dim, int spacedim>
1196 {
1197 return (const_cast<ParticleHandler<dim, spacedim> *>(this))->begin();
1198 }
1199
1200
1201
1202 template <int dim, int spacedim>
1205 {
1206 return particle_iterator(particle_container_owned_begin(),
1207 *property_pool,
1208 0);
1209 }
1210
1211
1212
1213 template <int dim, int spacedim>
1216 {
1217 return (const_cast<ParticleHandler<dim, spacedim> *>(this))->end();
1218 }
1219
1220
1221
1222 template <int dim, int spacedim>
1225 {
1226 return particle_iterator(particle_container_owned_end(), *property_pool, 0);
1227 }
1228
1229
1230
1231 template <int dim, int spacedim>
1234 {
1235 return (const_cast<ParticleHandler<dim, spacedim> *>(this))->begin_ghost();
1236 }
1237
1238
1239
1240 template <int dim, int spacedim>
1243 {
1244 return particle_iterator(particle_container_ghost_begin(),
1245 *property_pool,
1246 0);
1247 }
1248
1249
1250
1251 template <int dim, int spacedim>
1254 {
1255 return (const_cast<ParticleHandler<dim, spacedim> *>(this))->end_ghost();
1256 }
1257
1258
1259
1260 template <int dim, int spacedim>
1263 {
1264 return particle_iterator(particle_container_ghost_end(), *property_pool, 0);
1265 }
1266
1267
1268
1269 template <int dim, int spacedim>
1272 {
1273 // We should always have at least the three anchor entries in the list of
1274 // particles
1275 Assert(!particles.empty(), ExcInternalError());
1276 typename particle_container::iterator begin =
1277 const_cast<particle_container &>(particles).begin();
1278 return ++begin;
1279 }
1280
1281
1282
1283 template <int dim, int spacedim>
1286 {
1287 Assert(!particles.empty(), ExcInternalError());
1288 return owned_particles_end;
1289 }
1290
1291
1292
1293 template <int dim, int spacedim>
1296 {
1297 Assert(!particles.empty(), ExcInternalError());
1298 typename particle_container::iterator begin = owned_particles_end;
1299 return ++begin;
1300 }
1301
1302
1303
1304 template <int dim, int spacedim>
1307 {
1308 Assert(!particles.empty(), ExcInternalError());
1309 typename particle_container::iterator end =
1310 const_cast<particle_container &>(particles).end();
1311 return --end;
1312 }
1313
1314
1315
1316 template <int dim, int spacedim>
1317 template <class Archive>
1318 inline void
1319 ParticleHandler<dim, spacedim>::serialize(Archive &ar, const unsigned int)
1320 {
1321 // Note that we do not serialize the particle data itself. Instead we
1322 // use the serialization functionality of the triangulation class, because
1323 // this guarantees that data is immediately shipped to new processes if
1324 // the domain is distributed differently after resuming from a checkpoint.
1325 ar //&particles
1326 &global_number_of_particles &global_max_particles_per_cell
1327 & next_free_particle_index;
1328 }
1329
1330
1331
1332 template <int dim, int spacedim>
1333 template <class VectorType>
1334 inline std::enable_if_t<
1335 std::is_convertible<VectorType *, Function<spacedim> *>::value == false>
1337 const VectorType &input_vector,
1338 const bool displace_particles)
1339 {
1340 AssertDimension(input_vector.size(),
1341 get_next_free_particle_index() * spacedim);
1342 for (auto &p : *this)
1343 {
1344 Point<spacedim> new_point(displace_particles ? p.get_location() :
1345 Point<spacedim>());
1346 const auto id = p.get_id();
1347 for (unsigned int i = 0; i < spacedim; ++i)
1348 new_point[i] += input_vector[id * spacedim + i];
1349 p.set_location(new_point);
1350 }
1351 sort_particles_into_subdomains_and_cells();
1352 }
1353
1354
1355
1356 template <int dim, int spacedim>
1357 template <class VectorType>
1358 inline void
1360 VectorType &output_vector,
1361 const bool add_to_output_vector)
1362 {
1363 AssertDimension(output_vector.size(),
1364 get_next_free_particle_index() * spacedim);
1365 for (const auto &p : *this)
1366 {
1367 auto point = p.get_location();
1368 const auto id = p.get_id();
1369 if (add_to_output_vector)
1370 for (unsigned int i = 0; i < spacedim; ++i)
1371 output_vector[id * spacedim + i] += point[i];
1372 else
1373 for (unsigned int i = 0; i < spacedim; ++i)
1374 output_vector[id * spacedim + i] = point[i];
1375 }
1376 if (add_to_output_vector)
1377 output_vector.compress(VectorOperation::add);
1378 else
1379 output_vector.compress(VectorOperation::insert);
1380 }
1381
1382} // namespace Particles
1383
1385
1386#endif
Abstract base class for mapping classes.
Definition mapping.h:317
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
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 register_load_callback_function(const bool serialization)
std::function< std::size_t()> size_callback
std::enable_if_t< std::is_convertible< VectorType *, Function< spacedim > * >::value==false > set_particle_positions(const VectorType &input_vector, const bool displace_particles=true)
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
void unpack_callback(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() 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
particle_iterator_range particles_in_cell(const typename Triangulation< dim, spacedim >::active_cell_iterator &cell)
SmartPointer< const Triangulation< dim, spacedim >, ParticleHandler< dim, spacedim > > triangulation
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
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
SmartPointer< const Mapping< dim, spacedim >, ParticleHandler< dim, spacedim > > mapping
std::unique_ptr< GridTools::Cache< dim, spacedim > > triangulation_cache
types::particle_index get_next_free_particle_index() const
std::vector< char > pack_callback(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const typename Triangulation< dim, spacedim >::CellStatus status) const
IndexSet locally_owned_particle_ids() const
Definition point.h:112
#define DEAL_II_DEPRECATED
Definition config.h:172
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcInternalError()
boost::signals2::signal< void(const typename Particles::ParticleIterator< dim, spacedim > &particle, const typename Triangulation< dim, spacedim >::active_cell_iterator &cell)> particle_lost
const ::Triangulation< dim, spacedim > & tria