Reference documentation for deal.II version 9.6.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_accessor.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_accessor_h
16#define dealii_particles_particle_accessor_h
17
18#include <deal.II/base/config.h>
19
21
22#include <deal.II/grid/tria.h>
24
27
29
30namespace Particles
31{
32 // Forward declarations
33#ifndef DOXYGEN
34 template <int, int>
35 class ParticleIterator;
36 template <int, int>
37 class ParticleHandler;
38#endif
39
43 template <int dim, int spacedim = dim>
45 {
46 public:
81 {
85 ParticlesInCell() = default;
86
97
101 std::vector<typename PropertyPool<dim, spacedim>::Handle> particles;
102
107 };
108
112 using particle_container = std::list<ParticlesInCell>;
113
117 void *
119
120
124 const void *
126
147 void
148 set_location(const Point<spacedim> &new_location);
149
155 const Point<spacedim> &
157
181
202 void
203 set_reference_location(const Point<dim> &new_reference_location);
204
208 const Point<dim> &
210
214 void
216
221 get_id() const;
222
243
248 bool
250
271 void
272 set_properties(const std::vector<double> &new_properties);
273
294 void
296
320 void
321 set_properties(const Tensor<1, dim> &new_properties);
322
330 {
331 // The implementation is up here inside the class declaration because
332 // NVCC (at least in 12.5 and 12.6) otherwise produce a compile error:
333 //
334 // error: no declaration matches ‘::ArrayView<__remove_cv(const
335 // double)> ::Particles::ParticleAccessor<dim,
336 // spacedim>::get_properties()’
337 //
338 // See https://github.com/dealii/dealii/issues/17148
340
342 }
343
344
352
365 void
367
373 std::size_t
375
381
387 template <class Archive>
388 void
389 save(Archive &ar, const unsigned int version) const;
390
400 template <class Archive>
401 void
402 load(Archive &ar, const unsigned int version);
403
404#ifdef DOXYGEN
410 template <class Archive>
411 void
412 serialize(Archive &archive, const unsigned int version);
413#else
414 // This macro defines the serialize() method that is compatible with
415 // the templated save() and load() method that have been implemented.
416 BOOST_SERIALIZATION_SPLIT_MEMBER()
417#endif
418
422 void
424
428 void
430
434 bool
436
440 bool
442
447 state() const;
448
449 private:
454
462 const typename particle_container::iterator particles_in_cell,
464 const unsigned int particle_index_within_cell);
465
472 get_handle() const;
473
479
484 typename particle_container::iterator particles_in_cell;
485
490
495
496 // Make ParticleIterator a friend to allow it constructing
497 // ParticleAccessors.
498 template <int, int>
499 friend class ParticleIterator;
500 template <int, int>
501 friend class ParticleHandler;
502 };
503
504
505
506 template <int dim, int spacedim>
507 template <class Archive>
508 inline void
509 ParticleAccessor<dim, spacedim>::load(Archive &ar, const unsigned int)
510 {
511 unsigned int n_properties = 0;
512
513 Point<spacedim> location;
514 Point<dim> reference_location;
516 ar &location &reference_location &id &n_properties;
517
518 set_location(location);
519 set_reference_location(reference_location);
520 set_id(id);
521
522 if (n_properties > 0)
523 {
524 ArrayView<double> properties(get_properties());
525 Assert(
526 properties.size() == n_properties,
528 "This particle was serialized with " +
529 std::to_string(n_properties) +
530 " properties, but the new property handler provides space for " +
531 std::to_string(properties.size()) +
532 " properties. Deserializing a particle only works for matching property sizes."));
533
534 ar &boost::serialization::make_array(properties.data(), n_properties);
535 }
536 }
537
538
539
540 template <int dim, int spacedim>
541 template <class Archive>
542 inline void
543 ParticleAccessor<dim, spacedim>::save(Archive &ar, const unsigned int) const
544 {
545 unsigned int n_properties = 0;
546 if ((property_pool != nullptr) &&
548 n_properties = get_properties().size();
549
550 Point<spacedim> location = get_location();
551 Point<dim> reference_location = get_reference_location();
552 types::particle_index id = get_id();
553
554 ar &location &reference_location &id &n_properties;
555
556 if (n_properties > 0)
557 ar &boost::serialization::make_array(get_properties().data(),
558 n_properties);
559 }
560
561
562 // ------------------------- inline functions ------------------------------
563
564 template <int dim, int spacedim>
566 : particles_in_cell(typename particle_container::iterator())
567 , property_pool(nullptr)
568 , particle_index_within_cell(numbers::invalid_unsigned_int)
569 {}
570
571
572
573 template <int dim, int spacedim>
575 const typename particle_container::iterator particles_in_cell,
576 const PropertyPool<dim, spacedim> &property_pool,
577 const unsigned int particle_index_within_cell)
578 : particles_in_cell(particles_in_cell)
579 , property_pool(const_cast<PropertyPool<dim, spacedim> *>(&property_pool))
580 , particle_index_within_cell(particle_index_within_cell)
581 {}
582
583
584
585 template <int dim, int spacedim>
586 inline const void *
588 const void *data)
589 {
591
592 const types::particle_index *id_data =
593 static_cast<const types::particle_index *>(data);
594 set_id(*id_data++);
595 const double *pdata = reinterpret_cast<const double *>(id_data);
596
597 Point<spacedim> location;
598 for (unsigned int i = 0; i < spacedim; ++i)
599 location[i] = *pdata++;
600 set_location(location);
601
602 Point<dim> reference_location;
603 for (unsigned int i = 0; i < dim; ++i)
604 reference_location[i] = *pdata++;
605 set_reference_location(reference_location);
606
607 // See if there are properties to load
608 if (has_properties())
609 {
610 const ArrayView<double> particle_properties =
611 property_pool->get_properties(get_handle());
612 const unsigned int size = particle_properties.size();
613 for (unsigned int i = 0; i < size; ++i)
614 particle_properties[i] = *pdata++;
615 }
616
617 return static_cast<const void *>(pdata);
618 }
619
620
621
622 template <int dim, int spacedim>
623 inline void *
625 void *data) const
626 {
628
629 types::particle_index *id_data = static_cast<types::particle_index *>(data);
630 *id_data = get_id();
631 ++id_data;
632 double *pdata = reinterpret_cast<double *>(id_data);
633
634 // Write location
635 for (unsigned int i = 0; i < spacedim; ++i, ++pdata)
636 *pdata = get_location()[i];
637
638 // Write reference location
639 for (unsigned int i = 0; i < dim; ++i, ++pdata)
640 *pdata = get_reference_location()[i];
641
642 // Write properties
643 if (has_properties())
644 {
645 const ArrayView<double> particle_properties =
646 property_pool->get_properties(get_handle());
647 for (unsigned int i = 0; i < particle_properties.size(); ++i, ++pdata)
648 *pdata = particle_properties[i];
649 }
650
651 return static_cast<void *>(pdata);
652 }
653
654
655
656 template <int dim, int spacedim>
657 inline void
659 {
661
662 property_pool->set_location(get_handle(), new_loc);
663 }
664
665
666
667 template <int dim, int spacedim>
668 inline const Point<spacedim> &
670 {
672
673 return property_pool->get_location(get_handle());
674 }
675
676
677
678 template <int dim, int spacedim>
679 inline Point<spacedim> &
681 {
683
684 return property_pool->get_location(get_handle());
685 }
686
687
688
689 template <int dim, int spacedim>
690 inline void
692 const Point<dim> &new_loc)
693 {
695
696 property_pool->set_reference_location(get_handle(), new_loc);
697 }
698
699
700
701 template <int dim, int spacedim>
702 inline const Point<dim> &
704 {
706
707 return property_pool->get_reference_location(get_handle());
708 }
709
710
711
712 template <int dim, int spacedim>
715 {
717
718 return property_pool->get_id(get_handle());
719 }
720
721
722
723 template <int dim, int spacedim>
726 {
728
729 return get_handle();
730 }
731
732
733
734 template <int dim, int spacedim>
735 inline void
737 {
739
740 property_pool->set_id(get_handle(), new_id);
741 }
742
743
744
745 template <int dim, int spacedim>
746 inline bool
748 {
750
751 // Particles always have a property pool associated with them,
752 // but we can access properties only if there is a valid handle.
753 // The only way a particle can have no valid handle if it has
754 // been moved-from -- but that leaves an object in an invalid
755 // state, and so we can just assert that that can't be the case.
758 return (property_pool->n_properties_per_slot() > 0);
759 }
760
761
762
763 template <int dim, int spacedim>
764 inline void
766 const std::vector<double> &new_properties)
767 {
769
770 set_properties(
771 ArrayView<const double>(new_properties.data(), new_properties.size()));
772 }
773
774
775
776 template <int dim, int spacedim>
777 inline void
779 const ArrayView<const double> &new_properties)
780 {
782
783 const ArrayView<double> property_values =
784 property_pool->get_properties(get_handle());
785
786 Assert(new_properties.size() == property_values.size(),
788 "You are trying to assign properties with an incompatible length. "
789 "The particle has space to store " +
790 std::to_string(property_values.size()) +
791 " properties, but you are trying to assign " +
792 std::to_string(new_properties.size()) +
793 " properties. This is not allowed."));
794
795 if (property_values.size() > 0)
796 std::copy(new_properties.begin(),
797 new_properties.end(),
798 property_values.begin());
799 }
800
801
802
803 template <int dim, int spacedim>
804 inline void
806 const Tensor<1, dim> &new_properties)
807 {
809
810 // A Tensor object is not an array, so we cannot just create an
811 // ArrayView object for it. Rather, copy the data into a true
812 // array and make the ArrayView from that.
813 double array[dim];
814 for (unsigned int d = 0; d < dim; ++d)
815 array[d] = new_properties[d];
816
817 set_properties(make_array_view(array));
818 }
819
820
821
822 template <int dim, int spacedim>
825 {
827
828 return property_pool->get_properties(get_handle());
829 }
830
831
832
833 template <int dim, int spacedim>
834 inline void
836 PropertyPool<dim, spacedim> &new_property_pool)
837 {
838 Assert(&new_property_pool == property_pool, ExcInternalError());
839 (void)new_property_pool;
840 }
841
842
843
844 template <int dim, int spacedim>
845 inline const typename Triangulation<dim, spacedim>::cell_iterator &
847 {
849 Assert(particles_in_cell->cell.state() == IteratorState::valid,
851
852 return particles_in_cell->cell;
853 }
854
855
856
857 // template <int dim, int spacedim>
858 // inline ArrayView<double>
859 // ParticleAccessor<dim, spacedim>::get_properties()
860
861
862
863 template <int dim, int spacedim>
864 inline std::size_t
866 {
868
869 std::size_t size = sizeof(get_id()) +
870 sizeof(double) * spacedim + // get_location()
871 sizeof(double) * dim; // get_reference_location()
872
873 if (has_properties())
874 {
875 size += sizeof(double) * get_properties().size();
876 }
877 return size;
878 }
879
880
881
882 template <int dim, int spacedim>
883 inline void
885 {
887
888 ++particle_index_within_cell;
889
890 if (particle_index_within_cell >= particles_in_cell->particles.size())
891 {
892 particle_index_within_cell = 0;
893 ++particles_in_cell;
894 }
895 }
896
897
898
899 template <int dim, int spacedim>
900 inline void
902 {
904
905 if (particle_index_within_cell > 0)
906 --particle_index_within_cell;
907 else
908 {
909 --particles_in_cell;
910 particle_index_within_cell = particles_in_cell->particles.empty() ?
911 0 :
912 particles_in_cell->particles.size() - 1;
913 }
914 }
915
916
917
918 template <int dim, int spacedim>
919 inline bool
921 const ParticleAccessor<dim, spacedim> &other) const
922 {
923 return !(*this == other);
924 }
925
926
927
928 template <int dim, int spacedim>
929 inline bool
931 const ParticleAccessor<dim, spacedim> &other) const
932 {
933 return (property_pool == other.property_pool) &&
934 (particles_in_cell == other.particles_in_cell) &&
935 (particle_index_within_cell == other.particle_index_within_cell);
936 }
937
938
939
940 template <int dim, int spacedim>
943 {
944 if (property_pool != nullptr &&
945 particles_in_cell->cell.state() == IteratorState::valid &&
946 particle_index_within_cell < particles_in_cell->particles.size())
948 else if (property_pool != nullptr &&
949 particles_in_cell->cell.state() == IteratorState::past_the_end &&
950 particle_index_within_cell == 0)
952 else
954
956 }
957
958
959
960 template <int dim, int spacedim>
963 {
964 return particles_in_cell->particles[particle_index_within_cell];
965 }
966
967
968
969 template <int dim, int spacedim>
970 inline const typename PropertyPool<dim, spacedim>::Handle &
972 {
973 return particles_in_cell->particles[particle_index_within_cell];
974 }
975
976} // namespace Particles
977
979
980namespace boost
981{
982 namespace geometry
983 {
984 namespace index
985 {
986 // Forward declaration of bgi::indexable
987 template <class T>
988 struct indexable;
989
994 template <int dim, int spacedim>
995 struct indexable<::Particles::ParticleAccessor<dim, spacedim>>
996 {
1001 using result_type = const ::Point<spacedim> &;
1002
1004 operator()(const ::Particles::ParticleAccessor<dim, spacedim>
1005 &accessor) const
1006 {
1007 return accessor.get_location();
1008 }
1009 };
1010 } // namespace index
1011 } // namespace geometry
1012} // namespace boost
1013
1014#endif
ArrayView< std::remove_reference_t< typename std::iterator_traits< Iterator >::reference >, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
Definition array_view.h:949
iterator begin() const
Definition array_view.h:702
iterator end() const
Definition array_view.h:711
value_type * data() const noexcept
Definition array_view.h:661
std::size_t size() const
Definition array_view.h:684
particle_container::iterator particles_in_cell
const Point< spacedim > & get_location() const
types::particle_index get_id() const
void save(Archive &ar, const unsigned int version) const
types::particle_index get_local_index() const
std::list< ParticlesInCell > particle_container
const void * read_particle_data_from_memory(const void *data)
void load(Archive &ar, const unsigned int version)
const Triangulation< dim, spacedim >::cell_iterator & get_surrounding_cell() const
IteratorState::IteratorStates state() const
bool operator!=(const ParticleAccessor< dim, spacedim > &other) const
void set_properties(const Tensor< 1, dim > &new_properties)
void set_location(const Point< spacedim > &new_location)
PropertyPool< dim, spacedim >::Handle & get_handle()
void serialize(Archive &archive, const unsigned int version)
void set_property_pool(PropertyPool< dim, spacedim > &property_pool)
ArrayView< double > get_properties()
void set_reference_location(const Point< dim > &new_reference_location)
const Point< dim > & get_reference_location() const
std::size_t serialized_size_in_bytes() const
const PropertyPool< dim, spacedim >::Handle & get_handle() const
void set_properties(const std::vector< double > &new_properties)
Point< spacedim > & get_location()
ParticleAccessor(const typename particle_container::iterator particles_in_cell, const PropertyPool< dim, spacedim > &property_pool, const unsigned int particle_index_within_cell)
PropertyPool< dim, spacedim > * property_pool
void * write_particle_data_to_memory(void *data) const
void set_properties(const ArrayView< const double > &new_properties)
ArrayView< const double > get_properties() const
bool operator==(const ParticleAccessor< dim, spacedim > &other) const
void set_id(const types::particle_index &new_id)
ArrayView< double, ::MemorySpace::Host > get_properties(const Handle handle)
Definition point.h:111
#define DEAL_II_DEPRECATED
Definition config.h:207
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:503
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:504
#define Assert(cond, exc)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
@ past_the_end
Iterator reached end of container.
@ valid
Iterator points to a valid object.
@ invalid
Iterator is invalid, probably due to an error.
Triangulation< dim, spacedim >::active_cell_iterator cell
ParticlesInCell(const std::vector< typename PropertyPool< dim, spacedim >::Handle > &particles, const typename Triangulation< dim, spacedim >::active_cell_iterator &cell)
std::vector< typename PropertyPool< dim, spacedim >::Handle > particles
result_type operator()(const ::Particles::ParticleAccessor< dim, spacedim > &accessor) const