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_accessor.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2017 - 2022 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_accessor_h
17#define dealii_particles_particle_accessor_h
18
19#include <deal.II/base/config.h>
20
22
23#include <deal.II/grid/tria.h>
25
28
30
31namespace Particles
32{
33 // Forward declarations
34#ifndef DOXYGEN
35 template <int, int>
36 class ParticleIterator;
37 template <int, int>
38 class ParticleHandler;
39#endif
40
44 template <int dim, int spacedim = dim>
46 {
47 public:
82 {
86 ParticlesInCell() = default;
87
92 const std::vector<typename PropertyPool<dim, spacedim>::Handle>
93 &particles,
96 , cell(cell)
97 {}
98
102 std::vector<typename PropertyPool<dim, spacedim>::Handle> particles;
103
108 };
109
113 using particle_container = std::list<ParticlesInCell>;
114
118 void *
120
121
125 const void *
127
148 void
149 set_location(const Point<spacedim> &new_location);
150
156 const Point<spacedim> &
158
179 void
180 set_reference_location(const Point<dim> &new_reference_location);
181
185 const Point<dim> &
187
191 void
193
198 get_id() const;
199
220
225 bool
227
248 void
249 set_properties(const std::vector<double> &new_properties);
250
271 void
273
297 void
298 set_properties(const Tensor<1, dim> &new_properties);
299
307
315
328 void
330
336 std::size_t
338
347
356
362 template <class Archive>
363 void
364 save(Archive &ar, const unsigned int version) const;
365
375 template <class Archive>
376 void
377 load(Archive &ar, const unsigned int version);
378
379#ifdef DOXYGEN
385 template <class Archive>
386 void
387 serialize(Archive &archive, const unsigned int version);
388#else
389 // This macro defines the serialize() method that is compatible with
390 // the templated save() and load() method that have been implemented.
391 BOOST_SERIALIZATION_SPLIT_MEMBER()
392#endif
393
397 void
399
403 void
405
409 bool
411
415 bool
417
422 state() const;
423
424 private:
429
437 const typename particle_container::iterator particles_in_cell,
439 const unsigned int particle_index_within_cell);
440
447 get_handle() const;
448
454
459 typename particle_container::iterator particles_in_cell;
460
465
470
471 // Make ParticleIterator a friend to allow it constructing
472 // ParticleAccessors.
473 template <int, int>
474 friend class ParticleIterator;
475 template <int, int>
476 friend class ParticleHandler;
477 };
478
479
480
481 template <int dim, int spacedim>
482 template <class Archive>
483 inline void
484 ParticleAccessor<dim, spacedim>::load(Archive &ar, const unsigned int)
485 {
486 unsigned int n_properties = 0;
487
488 Point<spacedim> location;
489 Point<dim> reference_location;
491 ar &location &reference_location &id &n_properties;
492
493 set_location(location);
494 set_reference_location(reference_location);
495 set_id(id);
496
497 if (n_properties > 0)
498 {
499 ArrayView<double> properties(get_properties());
500 Assert(
501 properties.size() == n_properties,
503 "This particle was serialized with " +
504 std::to_string(n_properties) +
505 " properties, but the new property handler provides space for " +
506 std::to_string(properties.size()) +
507 " properties. Deserializing a particle only works for matching property sizes."));
508
509 ar &boost::serialization::make_array(properties.data(), n_properties);
510 }
511 }
512
513
514
515 template <int dim, int spacedim>
516 template <class Archive>
517 inline void
518 ParticleAccessor<dim, spacedim>::save(Archive &ar, const unsigned int) const
519 {
520 unsigned int n_properties = 0;
521 if ((property_pool != nullptr) &&
523 n_properties = get_properties().size();
524
525 Point<spacedim> location = get_location();
526 Point<dim> reference_location = get_reference_location();
527 types::particle_index id = get_id();
528
529 ar &location &reference_location &id &n_properties;
530
531 if (n_properties > 0)
532 ar &boost::serialization::make_array(get_properties().data(),
533 n_properties);
534 }
535
536
537 // ------------------------- inline functions ------------------------------
538
539 template <int dim, int spacedim>
541 : particles_in_cell(typename particle_container::iterator())
542 , property_pool(nullptr)
543 , particle_index_within_cell(numbers::invalid_unsigned_int)
544 {}
545
546
547
548 template <int dim, int spacedim>
550 const typename particle_container::iterator particles_in_cell,
551 const PropertyPool<dim, spacedim> & property_pool,
552 const unsigned int particle_index_within_cell)
553 : particles_in_cell(particles_in_cell)
554 , property_pool(const_cast<PropertyPool<dim, spacedim> *>(&property_pool))
555 , particle_index_within_cell(particle_index_within_cell)
556 {}
557
558
559
560 template <int dim, int spacedim>
561 inline const void *
563 const void *data)
564 {
566
567 const types::particle_index *id_data =
568 static_cast<const types::particle_index *>(data);
569 set_id(*id_data++);
570 const double *pdata = reinterpret_cast<const double *>(id_data);
571
572 Point<spacedim> location;
573 for (unsigned int i = 0; i < spacedim; ++i)
574 location(i) = *pdata++;
575 set_location(location);
576
577 Point<dim> reference_location;
578 for (unsigned int i = 0; i < dim; ++i)
579 reference_location(i) = *pdata++;
580 set_reference_location(reference_location);
581
582 // See if there are properties to load
583 if (has_properties())
584 {
585 const ArrayView<double> particle_properties =
586 property_pool->get_properties(get_handle());
587 const unsigned int size = particle_properties.size();
588 for (unsigned int i = 0; i < size; ++i)
589 particle_properties[i] = *pdata++;
590 }
591
592 return static_cast<const void *>(pdata);
593 }
594
595
596
597 template <int dim, int spacedim>
598 inline void *
600 void *data) const
601 {
603
604 types::particle_index *id_data = static_cast<types::particle_index *>(data);
605 *id_data = get_id();
606 ++id_data;
607 double *pdata = reinterpret_cast<double *>(id_data);
608
609 // Write location
610 for (unsigned int i = 0; i < spacedim; ++i, ++pdata)
611 *pdata = get_location()[i];
612
613 // Write reference location
614 for (unsigned int i = 0; i < dim; ++i, ++pdata)
615 *pdata = get_reference_location()[i];
616
617 // Write properties
618 if (has_properties())
619 {
620 const ArrayView<double> particle_properties =
621 property_pool->get_properties(get_handle());
622 for (unsigned int i = 0; i < particle_properties.size(); ++i, ++pdata)
623 *pdata = particle_properties[i];
624 }
625
626 return static_cast<void *>(pdata);
627 }
628
629
630
631 template <int dim, int spacedim>
632 inline void
634 {
636
637 property_pool->set_location(get_handle(), new_loc);
638 }
639
640
641
642 template <int dim, int spacedim>
643 inline const Point<spacedim> &
645 {
647
648 return property_pool->get_location(get_handle());
649 }
650
651
652
653 template <int dim, int spacedim>
654 inline void
656 const Point<dim> &new_loc)
657 {
659
660 property_pool->set_reference_location(get_handle(), new_loc);
661 }
662
663
664
665 template <int dim, int spacedim>
666 inline const Point<dim> &
668 {
670
671 return property_pool->get_reference_location(get_handle());
672 }
673
674
675
676 template <int dim, int spacedim>
679 {
681
682 return property_pool->get_id(get_handle());
683 }
684
685
686
687 template <int dim, int spacedim>
690 {
692
693 return get_handle();
694 }
695
696
697
698 template <int dim, int spacedim>
699 inline void
701 {
703
704 property_pool->set_id(get_handle(), new_id);
705 }
706
707
708
709 template <int dim, int spacedim>
710 inline bool
712 {
714
715 // Particles always have a property pool associated with them,
716 // but we can access properties only if there is a valid handle.
717 // The only way a particle can have no valid handle if it has
718 // been moved-from -- but that leaves an object in an invalid
719 // state, and so we can just assert that that can't be the case.
722 return (property_pool->n_properties_per_slot() > 0);
723 }
724
725
726
727 template <int dim, int spacedim>
728 inline void
730 const std::vector<double> &new_properties)
731 {
733
734 set_properties(
735 ArrayView<const double>(new_properties.data(), new_properties.size()));
736 }
737
738
739
740 template <int dim, int spacedim>
741 inline void
743 const ArrayView<const double> &new_properties)
744 {
746
747 const ArrayView<double> property_values =
748 property_pool->get_properties(get_handle());
749
750 Assert(new_properties.size() == property_values.size(),
752 "You are trying to assign properties with an incompatible length. "
753 "The particle has space to store " +
754 std::to_string(property_values.size()) +
755 " properties, but you are trying to assign " +
756 std::to_string(new_properties.size()) +
757 " properties. This is not allowed."));
758
759 if (property_values.size() > 0)
760 std::copy(new_properties.begin(),
761 new_properties.end(),
762 property_values.begin());
763 }
764
765
766
767 template <int dim, int spacedim>
768 inline void
770 const Tensor<1, dim> &new_properties)
771 {
773
774 // A Tensor object is not an array, so we cannot just create an
775 // ArrayView object for it. Rather, copy the data into a true
776 // array and make the ArrayView from that.
777 double array[dim];
778 for (unsigned int d = 0; d < dim; ++d)
779 array[d] = new_properties[d];
780
781 set_properties(make_array_view(array));
782 }
783
784
785
786 template <int dim, int spacedim>
787 inline const ArrayView<const double>
789 {
791
792 return property_pool->get_properties(get_handle());
793 }
794
795
796
797 template <int dim, int spacedim>
798 inline void
800 PropertyPool<dim, spacedim> &new_property_pool)
801 {
802 Assert(&new_property_pool == property_pool, ExcInternalError());
803 (void)new_property_pool;
804 }
805
806
807
808 template <int dim, int spacedim>
809 inline const typename Triangulation<dim, spacedim>::cell_iterator &
811 {
813 Assert(particles_in_cell->cell.state() == IteratorState::valid,
815
816 return particles_in_cell->cell;
817 }
818
819
820
821 template <int dim, int spacedim>
822 inline const typename Triangulation<dim, spacedim>::cell_iterator &
824 const Triangulation<dim, spacedim> & /*triangulation*/) const
825 {
827 Assert(particles_in_cell->cell.state() == IteratorState::valid,
829
830 return particles_in_cell->cell;
831 }
832
833
834
835 template <int dim, int spacedim>
836 inline const ArrayView<double>
838 {
840
841 return property_pool->get_properties(get_handle());
842 }
843
844
845
846 template <int dim, int spacedim>
847 inline std::size_t
849 {
851
852 std::size_t size = sizeof(get_id()) + sizeof(get_location()) +
853 sizeof(get_reference_location());
854
855 if (has_properties())
856 {
857 size += sizeof(double) * get_properties().size();
858 }
859 return size;
860 }
861
862
863
864 template <int dim, int spacedim>
865 inline void
867 {
869
870 ++particle_index_within_cell;
871
872 if (particle_index_within_cell >= particles_in_cell->particles.size())
873 {
874 particle_index_within_cell = 0;
875 ++particles_in_cell;
876 }
877 }
878
879
880
881 template <int dim, int spacedim>
882 inline void
884 {
886
887 if (particle_index_within_cell > 0)
888 --particle_index_within_cell;
889 else
890 {
891 --particles_in_cell;
892 particle_index_within_cell = particles_in_cell->particles.empty() ?
893 0 :
894 particles_in_cell->particles.size() - 1;
895 }
896 }
897
898
899
900 template <int dim, int spacedim>
901 inline bool
903 const ParticleAccessor<dim, spacedim> &other) const
904 {
905 return !(*this == other);
906 }
907
908
909
910 template <int dim, int spacedim>
911 inline bool
913 const ParticleAccessor<dim, spacedim> &other) const
914 {
915 return (property_pool == other.property_pool) &&
916 (particles_in_cell == other.particles_in_cell) &&
917 (particle_index_within_cell == other.particle_index_within_cell);
918 }
919
920
921
922 template <int dim, int spacedim>
925 {
926 if (property_pool != nullptr &&
927 particles_in_cell->cell.state() == IteratorState::valid &&
928 particle_index_within_cell < particles_in_cell->particles.size())
930 else if (property_pool != nullptr &&
931 particles_in_cell->cell.state() == IteratorState::past_the_end &&
932 particle_index_within_cell == 0)
934 else
936
938 }
939
940
941
942 template <int dim, int spacedim>
945 {
946 return particles_in_cell->particles[particle_index_within_cell];
947 }
948
949
950
951 template <int dim, int spacedim>
952 inline const typename PropertyPool<dim, spacedim>::Handle &
954 {
955 return particles_in_cell->particles[particle_index_within_cell];
956 }
957
958} // namespace Particles
959
961
962namespace boost
963{
964 namespace geometry
965 {
966 namespace index
967 {
968 // Forward declaration of bgi::indexable
969 template <class T>
970 struct indexable;
971
976 template <int dim, int spacedim>
977 struct indexable<::Particles::ParticleAccessor<dim, spacedim>>
978 {
983 using result_type = const ::Point<spacedim> &;
984
986 operator()(const ::Particles::ParticleAccessor<dim, spacedim>
987 &accessor) const
988 {
989 return accessor.get_location();
990 }
991 };
992 } // namespace index
993 } // namespace geometry
994} // namespace boost
995
996#endif
ArrayView< typename std::remove_reference< typename std::iterator_traits< Iterator >::reference >::type, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
Definition array_view.h:704
iterator begin() const
Definition array_view.h:594
iterator end() const
Definition array_view.h:603
value_type * data() const noexcept
Definition array_view.h:553
std::size_t size() const
Definition array_view.h:576
const Triangulation< dim, spacedim >::cell_iterator & get_surrounding_cell(const Triangulation< dim, spacedim > &triangulation) const
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)
const ArrayView< double > get_properties()
void set_location(const Point< spacedim > &new_location)
const ArrayView< const double > get_properties() const
PropertyPool< dim, spacedim >::Handle & get_handle()
void serialize(Archive &archive, const unsigned int version)
void set_property_pool(PropertyPool< dim, spacedim > &property_pool)
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)
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)
bool operator==(const ParticleAccessor< dim, spacedim > &other) const
void set_id(const types::particle_index &new_id)
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)
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.
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
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