Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-3087-ga35b476a78 2025-04-19 20:30:01+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
property_pool.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_property_pool_h
16#define dealii_particles_property_pool_h
17
18#include <deal.II/base/config.h>
19
21#include <deal.II/base/point.h>
22
23
25
26namespace types
27{
28 /* Type definitions */
29
30#ifdef DEAL_II_WITH_64BIT_INDICES
42 using particle_index = std::uint64_t;
43
44# ifdef DEAL_II_WITH_MPI
53# define DEAL_II_PARTICLE_INDEX_MPI_TYPE MPI_UINT64_T
54# endif
55
56#else
68 using particle_index = unsigned int;
69
70# ifdef DEAL_II_WITH_MPI
75# define DEAL_II_PARTICLE_INDEX_MPI_TYPE MPI_UNSIGNED
76# endif
77#endif
78} // namespace types
79
80namespace Particles
81{
104 template <int dim, int spacedim = dim>
106 {
107 public:
113 using Handle = unsigned int;
114
118 static const Handle invalid_handle;
119
123 PropertyPool(const unsigned int n_properties_per_slot);
124
131
138 void
139 clear();
140
146 Handle
148
153 void
155
160 const Point<spacedim> &
161 get_location(const Handle handle) const;
162
168 get_location(const Handle handle);
169
173 void
174 set_location(const Handle handle, const Point<spacedim> &new_location);
175
180 const Point<dim> &
181 get_reference_location(const Handle handle) const;
182
187 void
189 const Point<dim> &new_reference_location);
190
196 get_id(const Handle handle) const;
197
202 void
203 set_id(const Handle handle, const types::particle_index &new_id);
204
210 get_properties(const Handle handle)
211 {
212 // The implementation is up here inside the class declaration because
213 // NVCC (at least in 12.5 and 12.6) otherwise produce a compile error:
214 //
215 // error: no declaration matches ‘::ArrayView<__remove_cv(const
216 // double)> ::Particles::PropertyPool<dim,
217 // spacedim>::get_properties(Handle)’
218 //
219 // See https://github.com/dealii/dealii/issues/17148
220
221 const std::vector<double>::size_type data_index =
222 (handle != invalid_handle) ? handle * n_properties : 0;
223
224 // Ideally we would need to assert that 'handle' has not been deallocated
225 // by searching through 'currently_available_handles'. However, this
226 // is expensive and this function is performance critical, so instead
227 // just check against the array range, and rely on the fact
228 // that handles are invalidated when handed over to
229 // deallocate_properties_array().
230 Assert(data_index <= properties.size() - n_properties,
232 "Invalid property handle. This can happen if the "
233 "handle was duplicated and then one copy was deallocated "
234 "before trying to access the properties."));
235
236 return ArrayView<double>(properties.data() + data_index, n_properties);
237 }
238
239
244 void
245 reserve(const std::size_t size);
246
250 unsigned int
251 n_properties_per_slot() const;
252
257 unsigned int
258 n_slots() const;
259
263 unsigned int
264 n_registered_slots() const;
265
274 void
275 sort_memory_slots(const std::vector<Handle> &handles_to_sort);
276
277 private:
281 const unsigned int n_properties;
282
288 std::vector<Point<spacedim>> locations;
289
295 std::vector<Point<dim>> reference_locations;
296
302 std::vector<types::particle_index> ids;
303
309 std::vector<double> properties;
310
318 std::vector<Handle> currently_available_handles;
319 };
320
321
322
323 /* ---------------------- inline and template functions ------------------ */
324
325 template <int dim, int spacedim>
326 inline const Point<spacedim> &
328 {
329 const std::vector<double>::size_type data_index =
330 (handle != invalid_handle) ? handle : 0;
331
332 // Ideally we would need to assert that 'handle' has not been deallocated
333 // by searching through 'currently_available_handles'. However, this
334 // is expensive and this function is performance critical, so instead
335 // just check against the array range, and rely on the fact
336 // that handles are invalidated when handed over to
337 // deallocate_properties_array().
338 Assert(data_index <= locations.size() - 1,
339 ExcMessage("Invalid location handle. This can happen if the "
340 "handle was duplicated and then one copy was deallocated "
341 "before trying to access the properties."));
342
343 return locations[data_index];
344 }
345
346
347
348 template <int dim, int spacedim>
349 inline Point<spacedim> &
351 {
352 const std::vector<double>::size_type data_index =
353 (handle != invalid_handle) ? handle : 0;
354
355 // Ideally we would need to assert that 'handle' has not been deallocated
356 // by searching through 'currently_available_handles'. However, this
357 // is expensive and this function is performance critical, so instead
358 // just check against the array range, and rely on the fact
359 // that handles are invalidated when handed over to
360 // deallocate_properties_array().
361 Assert(data_index <= locations.size() - 1,
362 ExcMessage("Invalid location handle. This can happen if the "
363 "handle was duplicated and then one copy was deallocated "
364 "before trying to access the properties."));
365
366 return locations[data_index];
367 }
368
369
370
371 template <int dim, int spacedim>
372 inline void
374 const Point<spacedim> &new_location)
375 {
376 const std::vector<double>::size_type data_index =
377 (handle != invalid_handle) ? handle : 0;
378
379 // Ideally we would need to assert that 'handle' has not been deallocated
380 // by searching through 'currently_available_handles'. However, this
381 // is expensive and this function is performance critical, so instead
382 // just check against the array range, and rely on the fact
383 // that handles are invalidated when handed over to
384 // deallocate_properties_array().
385 Assert(data_index <= locations.size() - 1,
386 ExcMessage("Invalid location handle. This can happen if the "
387 "handle was duplicated and then one copy was deallocated "
388 "before trying to access the properties."));
389
390 locations[data_index] = new_location;
391 }
392
393
394
395 template <int dim, int spacedim>
396 inline const Point<dim> &
398 {
399 const std::vector<double>::size_type data_index =
400 (handle != invalid_handle) ? handle : 0;
401
402 // Ideally we would need to assert that 'handle' has not been deallocated
403 // by searching through 'currently_available_handles'. However, this
404 // is expensive and this function is performance critical, so instead
405 // just check against the array range, and rely on the fact
406 // that handles are invalidated when handed over to
407 // deallocate_properties_array().
408 Assert(data_index <= reference_locations.size() - 1,
409 ExcMessage("Invalid location handle. This can happen if the "
410 "handle was duplicated and then one copy was deallocated "
411 "before trying to access the properties."));
412
413 return reference_locations[data_index];
414 }
415
416
417
418 template <int dim, int spacedim>
419 inline void
421 const Handle handle,
422 const Point<dim> &new_reference_location)
423 {
424 const std::vector<double>::size_type data_index =
425 (handle != invalid_handle) ? handle : 0;
426
427 // Ideally we would need to assert that 'handle' has not been deallocated
428 // by searching through 'currently_available_handles'. However, this
429 // is expensive and this function is performance critical, so instead
430 // just check against the array range, and rely on the fact
431 // that handles are invalidated when handed over to
432 // deallocate_properties_array().
433 Assert(data_index <= reference_locations.size() - 1,
434 ExcMessage("Invalid location handle. This can happen if the "
435 "handle was duplicated and then one copy was deallocated "
436 "before trying to access the properties."));
437
438 reference_locations[data_index] = new_reference_location;
439 }
440
441
442
443 template <int dim, int spacedim>
446 {
447 const std::vector<double>::size_type data_index =
448 (handle != invalid_handle) ? handle : 0;
449
450 // Ideally we would need to assert that 'handle' has not been deallocated
451 // by searching through 'currently_available_handles'. However, this
452 // is expensive and this function is performance critical, so instead
453 // just check against the array range, and rely on the fact
454 // that handles are invalidated when handed over to
455 // deallocate_properties_array().
456 Assert(data_index <= ids.size() - 1,
457 ExcMessage("Invalid location handle. This can happen if the "
458 "handle was duplicated and then one copy was deallocated "
459 "before trying to access the properties."));
460
461 return ids[data_index];
462 }
463
464
465
466 template <int dim, int spacedim>
467 inline void
469 const types::particle_index &new_id)
470 {
471 const std::vector<double>::size_type data_index =
472 (handle != invalid_handle) ? handle : 0;
473
474 // Ideally we would need to assert that 'handle' has not been deallocated
475 // by searching through 'currently_available_handles'. However, this
476 // is expensive and this function is performance critical, so instead
477 // just check against the array range, and rely on the fact
478 // that handles are invalidated when handed over to
479 // deallocate_properties_array().
480 Assert(data_index <= ids.size() - 1,
481 ExcMessage("Invalid location handle. This can happen if the "
482 "handle was duplicated and then one copy was deallocated "
483 "before trying to access the properties."));
484
485 ids[data_index] = new_id;
486 }
487
488
489
490 // template <int dim, int spacedim>
491 // inline ArrayView<double, ::MemorySpace::Host>
492 // PropertyPool<dim, spacedim>::get_properties(const Handle handle)
493
494
495
496 template <int dim, int spacedim>
497 inline unsigned int
499 {
500 return locations.size();
501 }
502
503
504} // namespace Particles
505
507
508#endif
const Point< spacedim > & get_location(const Handle handle) const
std::vector< Point< dim > > reference_locations
unsigned int n_properties_per_slot() const
unsigned int n_registered_slots() const
unsigned int n_slots() const
const unsigned int n_properties
void deregister_particle(Handle &handle)
void set_id(const Handle handle, const types::particle_index &new_id)
types::particle_index get_id(const Handle handle) const
void set_reference_location(const Handle handle, const Point< dim > &new_reference_location)
void reserve(const std::size_t size)
void set_location(const Handle handle, const Point< spacedim > &new_location)
const Point< dim > & get_reference_location(const Handle handle) const
void sort_memory_slots(const std::vector< Handle > &handles_to_sort)
Point< spacedim > & get_location(const Handle handle)
std::vector< Point< spacedim > > locations
std::vector< types::particle_index > ids
ArrayView< double, ::MemorySpace::Host > get_properties(const Handle handle)
static const Handle invalid_handle
std::vector< Handle > currently_available_handles
std::vector< double > properties
Definition point.h:113
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:35
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
#define Assert(cond, exc)
static ::ExceptionBase & ExcMessage(std::string arg1)
std::size_t size
Definition mpi.cc:745
Definition types.h:32
unsigned int particle_index