From: Rene Gassmoeller Date: Wed, 13 Sep 2017 15:18:10 +0000 (+0200) Subject: Add accessor and iterator classes X-Git-Tag: v9.0.0-rc1~1013^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=d133835cc0258a8772b66c307ae63949bad23e3f;p=dealii.git Add accessor and iterator classes --- diff --git a/include/deal.II/particles/particle.h b/include/deal.II/particles/particle.h index 9a2de57d11..d27082cb92 100644 --- a/include/deal.II/particles/particle.h +++ b/include/deal.II/particles/particle.h @@ -13,8 +13,8 @@ // // --------------------------------------------------------------------- -#ifndef dealii__particles_particle_h -#define dealii__particles_particle_h +#ifndef dealii_particles_particle_h +#define dealii_particles_particle_h #include @@ -37,6 +37,12 @@ namespace Particles { /* Type definitions */ + /** + * Typedef of cell level/index pair. TODO: replace this by the + * active_cell_index from deal.II 8.3 onwards. + */ + typedef std::pair LevelInd; + #ifdef DEAL_II_WITH_64BIT_INDICES /** * The type used for indices of particles. While in diff --git a/include/deal.II/particles/particle_accessor.h b/include/deal.II/particles/particle_accessor.h new file mode 100644 index 0000000000..94c51422c5 --- /dev/null +++ b/include/deal.II/particles/particle_accessor.h @@ -0,0 +1,205 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#ifndef dealii_particles_particle_accessor_h +#define dealii_particles_particle_accessor_h + +#include +#include +#include + +DEAL_II_NAMESPACE_OPEN + +namespace Particles +{ + template class ParticleIterator; + template class ParticleHandler; + + template + class ParticleAccessor + { + public: + /** + * Write particle data into a data array. The array is expected + * to be large enough to take the data, and the void pointer should + * point to the first element in which the data should be written. This + * function is meant for serializing all particle properties and + * afterwards de-serializing the properties by calling the appropriate + * constructor Particle(void *&data, PropertyPool *property_pool = NULL); + * + * @param [in,out] data The memory location to write particle data + * into. This pointer points to the begin of the memory, in which the + * data will be written and it will be advanced by the serialized size + * of this particle. + */ + void + write_data(void *&data) const; + + /** + * Set the location of this particle. Note that this does not check + * whether this is a valid location in the simulation domain. + * + * @param [in] new_location The new location for this particle. + */ + void + set_location (const Point &new_location); + + /** + * Get the location of this particle. + * + * @return The location of this particle. + */ + const Point & + get_location () const; + + /** + * Set the reference location of this particle. + * + * @param [in] new_reference_location The new reference location for + * this particle. + */ + void + set_reference_location (const Point &new_reference_location); + + /** + * Return the reference location of this particle in its current cell. + */ + const Point & + get_reference_location () const; + + /** + * Return the ID number of this particle. + */ + types::particle_index + get_id () const; + + /** + * Tell the particle where to store its properties (even if it does not + * own properties). Usually this is only done once per particle, but + * since the particle generator does not know about the properties + * we want to do it not at construction time. Another use for this + * function is after particle transfer to a new process. + */ + void + set_property_pool(PropertyPool &property_pool); + + /** + * Returns whether this particle has a valid property pool and a valid + * handle to properties. + */ + bool + has_properties () const; + + /** + * Set the properties of this particle. + * + * @param [in] new_properties A vector containing the + * new properties for this particle. + */ + void + set_properties (const std::vector &new_properties); + + /** + * Get write-access to properties of this particle. + * + * @return An ArrayView of the properties of this particle. + */ + const ArrayView + get_properties (); + + /** + * Get read-access to properties of this particle. + * + * @return An ArrayView of the properties of this particle. + */ + const ArrayView + get_properties () const; + + /** + * Returns the size in bytes this particle occupies if all of its data is + * serialized (i.e. the number of bytes that is written by the write_data + * function of this class). + */ + std::size_t + size() const; + + /** + * Get a cell iterator to the cell surrounding the current particle. + * As particles are organized in the structure of a triangulation, + * but the triangulation itself is not stored in the particle this + * operation requires a reference to the triangulation. + */ + typename parallel::distributed::Triangulation::cell_iterator + get_surrounding_cell (const parallel::distributed::Triangulation &triangulation) const; + + /** + * Serialize the contents of this class. + */ + template + void serialize (Archive &ar, const unsigned int version); + + /** + * Advance the ParticleAccessor to the next particle. + */ + void next (); + + /** + * Move the ParticleAccessor to the previous particle. + */ + void prev (); + + /** + * Inequality operator. + */ + bool operator != (const ParticleAccessor &other) const; + + /** + * Equality operator. + */ + bool operator == (const ParticleAccessor &other) const; + + protected: + /** + * Construct an accessor from a reference to a map and an iterator to the map. + * This constructor is protected so that it can only be accessed by friend + * classes. + */ + ParticleAccessor (const std::multimap > &map, + const typename std::multimap >::iterator &particle); + + private: + /** + * A pointer to the container that stores the particles. Obviously, + * this accessor is invalidated if the container changes. + */ + std::multimap > *map; + + /** + * An iterator into the container of particles. Obviously, + * this accessor is invalidated if the container changes. + */ + typename std::multimap >::iterator particle; + + /** + * Make ParticleIterator a friend to allow it constructing ParticleAccessors. + */ + template friend class ParticleIterator; + template friend class ParticleHandler; + }; +} + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/include/deal.II/particles/particle_iterator.h b/include/deal.II/particles/particle_iterator.h new file mode 100644 index 0000000000..83941dd6b6 --- /dev/null +++ b/include/deal.II/particles/particle_iterator.h @@ -0,0 +1,113 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#ifndef dealii_particles_particle_iterator_h +#define dealii_particles_particle_iterator_h + +#include + +DEAL_II_NAMESPACE_OPEN + +namespace Particles +{ + /** + * A class that is used to iterate over particles. Together with the + * ParticleAccessor class this is used to hide the internal implementation + * of the particle class and the particle container. + */ + template + class ParticleIterator: public std::iterator > + { + public: + /** + * Constructor of the iterator. Takes a reference to the particle + * container, and an iterator the the cell-particle pair. + */ + ParticleIterator (const std::multimap > &map, + const typename std::multimap >::iterator &particle); + + /** + * Dereferencing operator, returns a reference to an accessor. Usage is thus + * like (*i).get_id (); + */ + const ParticleAccessor &operator * () const; + + /** + * Dereferencing operator, non-@p const version. + */ + ParticleAccessor &operator * (); + + /** + * Dereferencing operator, returns a pointer of the particle pointed to. Usage + * is thus like i->get_id (); + * + * There is a @p const and a non-@p const version. + */ + const ParticleAccessor *operator -> () const; + + /** + * Dereferencing operator, non-@p const version. + */ + ParticleAccessor *operator -> (); + + /** + * Compare for equality. + */ + bool operator == (const ParticleIterator &) const; + + /** + * Compare for inequality. + */ + bool operator != (const ParticleIterator &) const; + + /** + * Prefix ++ operator: ++iterator. This operator advances + * the iterator to the next element and returns a reference to + * *this. + */ + ParticleIterator &operator ++ (); + + /** + * Postfix ++ operator: iterator++. This operator advances + * the iterator to the next element, but returns an iterator to the element + * previously pointed to. + */ + ParticleIterator operator ++ (int); + + /** + * Prefix -- operator: --iterator. This operator moves + * the iterator to the previous element and returns a reference to + * *this. + */ + ParticleIterator &operator -- (); + + /** + * Postfix -- operator: iterator--. This operator moves + * the iterator to the previous element, but returns an iterator to the element + * previously pointed to. + */ + ParticleIterator operator -- (int); + + private: + /** + * The accessor to the actual particle. + */ + ParticleAccessor accessor; + }; +} + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/source/particles/CMakeLists.txt b/source/particles/CMakeLists.txt index f355e5050f..8eaa91767d 100644 --- a/source/particles/CMakeLists.txt +++ b/source/particles/CMakeLists.txt @@ -18,11 +18,15 @@ INCLUDE_DIRECTORIES(BEFORE ${CMAKE_CURRENT_BINARY_DIR}) SET(_src particle.cc + particle_accessor.cc + particle_iterator.cc property_pool.cc ) SET(_inst particle.inst.in + particle_accessor.inst.in + particle_iterator.inst.in property_pool.inst.in ) diff --git a/source/particles/particle.cc b/source/particles/particle.cc index e2917c5f00..b898f21193 100644 --- a/source/particles/particle.cc +++ b/source/particles/particle.cc @@ -215,8 +215,8 @@ namespace Particles Particle::size () const { std::size_t size = sizeof(types::particle_index) - + sizeof(location) - + sizeof(reference_location); + + sizeof(location) + + sizeof(reference_location); if (has_properties()) { @@ -287,7 +287,7 @@ namespace Particles Particle::has_properties () const { return (property_pool != NULL) - && (properties != PropertyPool::invalid_handle); + && (properties != PropertyPool::invalid_handle); } @@ -297,8 +297,8 @@ namespace Particles Particle::set_properties (const ArrayView &new_properties) { Assert(property_pool != NULL, - ExcMessage("A particle was asked to set its properties without an " - "associated PropertyPool.")) + ExcMessage("A particle was asked to set its properties without an " + "associated PropertyPool.")) if (properties == PropertyPool::invalid_handle) properties = property_pool->allocate_properties_array(); @@ -316,7 +316,7 @@ namespace Particles { Assert(has_properties(), ExcMessage("A particle without additional properties was asked for " - "its properties.")); + "its properties.")); return property_pool->get_properties(properties); } @@ -329,7 +329,7 @@ namespace Particles { Assert(has_properties(), ExcMessage("A particle without additional properties was asked for " - "its properties.")); + "its properties.")); return property_pool->get_properties(properties); } diff --git a/source/particles/particle.inst.in b/source/particles/particle.inst.in index 9ea1b7927d..1f1327202a 100644 --- a/source/particles/particle.inst.in +++ b/source/particles/particle.inst.in @@ -16,9 +16,11 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS; number : REAL_SCALARS) { +#if deal_II_dimension <= deal_II_space_dimension namespace Particles \{ template class Particle ; \} +#endif } diff --git a/source/particles/particle_accessor.cc b/source/particles/particle_accessor.cc new file mode 100644 index 0000000000..16bbd14c65 --- /dev/null +++ b/source/particles/particle_accessor.cc @@ -0,0 +1,223 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#include + +DEAL_II_NAMESPACE_OPEN + +namespace Particles +{ + template + ParticleAccessor::ParticleAccessor (const std::multimap > &map, + const typename std::multimap >::iterator &particle) + : + map (const_cast > *> (&map)), + particle (particle) + {} + + + + template + void + ParticleAccessor::write_data (void *&data) const + { + Assert(particle != map->end(), + ExcInternalError()); + + particle->second.write_data(data); + } + + + + template + void + ParticleAccessor::set_location (const Point &new_loc) + { + Assert(particle != map->end(), + ExcInternalError()); + + particle->second.set_location(new_loc); + } + + + + template + const Point & + ParticleAccessor::get_location () const + { + Assert(particle != map->end(), + ExcInternalError()); + + return particle->second.get_location(); + } + + + + template + void + ParticleAccessor::set_reference_location (const Point &new_loc) + { + Assert(particle != map->end(), + ExcInternalError()); + + particle->second.set_reference_location(new_loc); + } + + + + template + const Point & + ParticleAccessor::get_reference_location () const + { + Assert(particle != map->end(), + ExcInternalError()); + + return particle->second.get_reference_location(); + } + + + + template + types::particle_index + ParticleAccessor::get_id () const + { + Assert(particle != map->end(), + ExcInternalError()); + + return particle->second.get_id(); + } + + + + template + void + ParticleAccessor::set_property_pool (PropertyPool &new_property_pool) + { + Assert(particle != map->end(), + ExcInternalError()); + + particle->second.set_property_pool(new_property_pool); + } + + + + template + void + ParticleAccessor::set_properties (const std::vector &new_properties) + { + Assert(particle != map->end(), + ExcInternalError()); + + particle->second.set_properties(ArrayView(&(new_properties[0]),new_properties.size())); + return; + } + + + + template + const ArrayView + ParticleAccessor::get_properties () const + { + Assert(particle != map->end(), + ExcInternalError()); + + return particle->second.get_properties(); + } + + + + template + typename parallel::distributed::Triangulation::cell_iterator + ParticleAccessor::get_surrounding_cell (const parallel::distributed::Triangulation &triangulation) const + { + Assert(particle != map->end(), + ExcInternalError()); + + const typename parallel::distributed::Triangulation::cell_iterator cell (&triangulation, + particle->first.first, + particle->first.second); + return cell; + } + + + + template + const ArrayView + ParticleAccessor::get_properties () + { + Assert(particle != map->end(), + ExcInternalError()); + + return particle->second.get_properties(); + } + + + + template + std::size_t + ParticleAccessor::size () const + { + Assert(particle != map->end(), + ExcInternalError()); + + return particle->second.size(); + } + + + + template + void + ParticleAccessor::next () + { + Assert (particle != map->end(),ExcInternalError()); + ++particle; + } + + + + template + void + ParticleAccessor::prev () + { + Assert (particle != map->begin(),ExcInternalError()); + --particle; + } + + + + template + bool + ParticleAccessor::operator != (const ParticleAccessor &other) const + { + return (map != other.map) || (particle != other.particle); + } + + + + template + bool + ParticleAccessor::operator == (const ParticleAccessor &other) const + { + return (map == other.map) && (particle == other.particle); + } +} + +DEAL_II_NAMESPACE_CLOSE + +DEAL_II_NAMESPACE_OPEN + +#include "particle_accessor.inst" + +DEAL_II_NAMESPACE_CLOSE diff --git a/source/particles/particle_accessor.inst.in b/source/particles/particle_accessor.inst.in new file mode 100644 index 0000000000..71a9286f7c --- /dev/null +++ b/source/particles/particle_accessor.inst.in @@ -0,0 +1,26 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) +{ +#if deal_II_dimension <= deal_II_space_dimension + namespace Particles + \{ + template + class ParticleAccessor ; + \} +#endif +} diff --git a/source/particles/particle_iterator.cc b/source/particles/particle_iterator.cc new file mode 100644 index 0000000000..4fa7873382 --- /dev/null +++ b/source/particles/particle_iterator.cc @@ -0,0 +1,135 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#include + +DEAL_II_NAMESPACE_OPEN + +namespace Particles +{ + template + ParticleIterator::ParticleIterator (const std::multimap > &map, + const typename std::multimap >::iterator &particle) + : + accessor (map, particle) + {} + + + + template + ParticleAccessor & + ParticleIterator::operator *() + { + return accessor; + } + + + + template + ParticleAccessor * + ParticleIterator::operator ->() + { + return &(this->operator* ()); + } + + + + template + const ParticleAccessor & + ParticleIterator::operator *() const + { + return accessor; + } + + + + template + const ParticleAccessor * + ParticleIterator::operator ->() const + { + return &(this->operator* ()); + } + + + + template + bool + ParticleIterator::operator != (const ParticleIterator &other) const + { + return accessor != other.accessor; + } + + + + template + bool + ParticleIterator::operator == (const ParticleIterator &other) const + { + return accessor == other.accessor; + } + + + + template + ParticleIterator & + ParticleIterator::operator++() + { + accessor.next(); + return *this; + } + + + + template + ParticleIterator + ParticleIterator::operator++(int) + { + ParticleIterator tmp(*this); + operator++ (); + + return tmp; + } + + + + template + ParticleIterator & + ParticleIterator::operator--() + { + accessor.prev(); + return *this; + } + + + + template + ParticleIterator + ParticleIterator::operator--(int) + { + ParticleIterator tmp(*this); + operator-- (); + + return tmp; + } +} + + +DEAL_II_NAMESPACE_CLOSE + +DEAL_II_NAMESPACE_OPEN + +#include "particle_iterator.inst" + +DEAL_II_NAMESPACE_CLOSE diff --git a/source/particles/particle_iterator.inst.in b/source/particles/particle_iterator.inst.in new file mode 100644 index 0000000000..2d6e7c9678 --- /dev/null +++ b/source/particles/particle_iterator.inst.in @@ -0,0 +1,26 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) +{ +#if deal_II_dimension <= deal_II_space_dimension + namespace Particles + \{ + template + class ParticleIterator ; + \} +#endif +}