//
// ---------------------------------------------------------------------
-#ifndef dealii__particles_particle_h
-#define dealii__particles_particle_h
+#ifndef dealii_particles_particle_h
+#define dealii_particles_particle_h
#include <deal.II/particles/property_pool.h>
{
/* 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<int, int> LevelInd;
+
#ifdef DEAL_II_WITH_64BIT_INDICES
/**
* The type used for indices of particles. While in
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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 <deal.II/particles/particle.h>
+#include <deal.II/base/array_view.h>
+#include <deal.II/distributed/tria.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace Particles
+{
+ template <int, int> class ParticleIterator;
+ template <int, int> class ParticleHandler;
+
+ template<int dim, int spacedim=dim>
+ 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<spacedim> &new_location);
+
+ /**
+ * Get the location of this particle.
+ *
+ * @return The location of this particle.
+ */
+ const Point<spacedim> &
+ 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<dim> &new_reference_location);
+
+ /**
+ * Return the reference location of this particle in its current cell.
+ */
+ const Point<dim> &
+ 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<double> &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<double> &new_properties);
+
+ /**
+ * Get write-access to properties of this particle.
+ *
+ * @return An ArrayView of the properties of this particle.
+ */
+ const ArrayView<double>
+ get_properties ();
+
+ /**
+ * Get read-access to properties of this particle.
+ *
+ * @return An ArrayView of the properties of this particle.
+ */
+ const ArrayView<const double>
+ 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<dim,spacedim>::cell_iterator
+ get_surrounding_cell (const parallel::distributed::Triangulation<dim,spacedim> &triangulation) const;
+
+ /**
+ * Serialize the contents of this class.
+ */
+ template <class Archive>
+ 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<dim,spacedim> &other) const;
+
+ /**
+ * Equality operator.
+ */
+ bool operator == (const ParticleAccessor<dim,spacedim> &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<types::LevelInd, Particle<dim,spacedim> > &map,
+ const typename std::multimap<types::LevelInd, Particle<dim,spacedim> >::iterator &particle);
+
+ private:
+ /**
+ * A pointer to the container that stores the particles. Obviously,
+ * this accessor is invalidated if the container changes.
+ */
+ std::multimap<types::LevelInd, Particle<dim,spacedim> > *map;
+
+ /**
+ * An iterator into the container of particles. Obviously,
+ * this accessor is invalidated if the container changes.
+ */
+ typename std::multimap<types::LevelInd, Particle<dim,spacedim> >::iterator particle;
+
+ /**
+ * Make ParticleIterator a friend to allow it constructing ParticleAccessors.
+ */
+ template <int, int> friend class ParticleIterator;
+ template <int, int> friend class ParticleHandler;
+ };
+}
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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/particles/particle_accessor.h>
+
+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<int dim, int spacedim=dim>
+ class ParticleIterator: public std::iterator<std::bidirectional_iterator_tag,ParticleAccessor<dim,spacedim> >
+ {
+ public:
+ /**
+ * Constructor of the iterator. Takes a reference to the particle
+ * container, and an iterator the the cell-particle pair.
+ */
+ ParticleIterator (const std::multimap<types::LevelInd, Particle<dim,spacedim> > &map,
+ const typename std::multimap<types::LevelInd, Particle<dim,spacedim> >::iterator &particle);
+
+ /**
+ * Dereferencing operator, returns a reference to an accessor. Usage is thus
+ * like <tt>(*i).get_id ();</tt>
+ */
+ const ParticleAccessor<dim,spacedim> &operator * () const;
+
+ /**
+ * Dereferencing operator, non-@p const version.
+ */
+ ParticleAccessor<dim,spacedim> &operator * ();
+
+ /**
+ * Dereferencing operator, returns a pointer of the particle pointed to. Usage
+ * is thus like <tt>i->get_id ();</tt>
+ *
+ * There is a @p const and a non-@p const version.
+ */
+ const ParticleAccessor<dim,spacedim> *operator -> () const;
+
+ /**
+ * Dereferencing operator, non-@p const version.
+ */
+ ParticleAccessor<dim,spacedim> *operator -> ();
+
+ /**
+ * Compare for equality.
+ */
+ bool operator == (const ParticleIterator<dim,spacedim> &) const;
+
+ /**
+ * Compare for inequality.
+ */
+ bool operator != (const ParticleIterator<dim,spacedim> &) const;
+
+ /**
+ * Prefix <tt>++</tt> operator: <tt>++iterator</tt>. This operator advances
+ * the iterator to the next element and returns a reference to
+ * <tt>*this</tt>.
+ */
+ ParticleIterator &operator ++ ();
+
+ /**
+ * Postfix <tt>++</tt> operator: <tt>iterator++</tt>. This operator advances
+ * the iterator to the next element, but returns an iterator to the element
+ * previously pointed to.
+ */
+ ParticleIterator operator ++ (int);
+
+ /**
+ * Prefix <tt>--</tt> operator: <tt>--iterator</tt>. This operator moves
+ * the iterator to the previous element and returns a reference to
+ * <tt>*this</tt>.
+ */
+ ParticleIterator &operator -- ();
+
+ /**
+ * Postfix <tt>--</tt> operator: <tt>iterator--</tt>. 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<dim,spacedim> accessor;
+ };
+}
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
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
)
Particle<dim,spacedim,PropertyType>::size () const
{
std::size_t size = sizeof(types::particle_index)
- + sizeof(location)
- + sizeof(reference_location);
+ + sizeof(location)
+ + sizeof(reference_location);
if (has_properties())
{
Particle<dim,spacedim,PropertyType>::has_properties () const
{
return (property_pool != NULL)
- && (properties != PropertyPool<PropertyType>::invalid_handle);
+ && (properties != PropertyPool<PropertyType>::invalid_handle);
}
Particle<dim,spacedim,PropertyType>::set_properties (const ArrayView<const PropertyType> &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<PropertyType>::invalid_handle)
properties = property_pool->allocate_properties_array();
{
Assert(has_properties(),
ExcMessage("A particle without additional properties was asked for "
- "its properties."));
+ "its properties."));
return property_pool->get_properties(properties);
}
{
Assert(has_properties(),
ExcMessage("A particle without additional properties was asked for "
- "its properties."));
+ "its properties."));
return property_pool->get_properties(properties);
}
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 <deal_II_dimension,deal_II_space_dimension,number>;
\}
+#endif
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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/particles/particle_accessor.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace Particles
+{
+ template <int dim, int spacedim>
+ ParticleAccessor<dim,spacedim>::ParticleAccessor (const std::multimap<types::LevelInd, Particle<dim,spacedim> > &map,
+ const typename std::multimap<types::LevelInd, Particle<dim,spacedim> >::iterator &particle)
+ :
+ map (const_cast<std::multimap<types::LevelInd, Particle<dim,spacedim> > *> (&map)),
+ particle (particle)
+ {}
+
+
+
+ template <int dim, int spacedim>
+ void
+ ParticleAccessor<dim,spacedim>::write_data (void *&data) const
+ {
+ Assert(particle != map->end(),
+ ExcInternalError());
+
+ particle->second.write_data(data);
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ ParticleAccessor<dim,spacedim>::set_location (const Point<spacedim> &new_loc)
+ {
+ Assert(particle != map->end(),
+ ExcInternalError());
+
+ particle->second.set_location(new_loc);
+ }
+
+
+
+ template <int dim, int spacedim>
+ const Point<spacedim> &
+ ParticleAccessor<dim,spacedim>::get_location () const
+ {
+ Assert(particle != map->end(),
+ ExcInternalError());
+
+ return particle->second.get_location();
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ ParticleAccessor<dim,spacedim>::set_reference_location (const Point<dim> &new_loc)
+ {
+ Assert(particle != map->end(),
+ ExcInternalError());
+
+ particle->second.set_reference_location(new_loc);
+ }
+
+
+
+ template <int dim, int spacedim>
+ const Point<dim> &
+ ParticleAccessor<dim,spacedim>::get_reference_location () const
+ {
+ Assert(particle != map->end(),
+ ExcInternalError());
+
+ return particle->second.get_reference_location();
+ }
+
+
+
+ template <int dim, int spacedim>
+ types::particle_index
+ ParticleAccessor<dim,spacedim>::get_id () const
+ {
+ Assert(particle != map->end(),
+ ExcInternalError());
+
+ return particle->second.get_id();
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ ParticleAccessor<dim,spacedim>::set_property_pool (PropertyPool<double> &new_property_pool)
+ {
+ Assert(particle != map->end(),
+ ExcInternalError());
+
+ particle->second.set_property_pool(new_property_pool);
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ ParticleAccessor<dim,spacedim>::set_properties (const std::vector<double> &new_properties)
+ {
+ Assert(particle != map->end(),
+ ExcInternalError());
+
+ particle->second.set_properties(ArrayView<const double>(&(new_properties[0]),new_properties.size()));
+ return;
+ }
+
+
+
+ template <int dim, int spacedim>
+ const ArrayView<const double>
+ ParticleAccessor<dim,spacedim>::get_properties () const
+ {
+ Assert(particle != map->end(),
+ ExcInternalError());
+
+ return particle->second.get_properties();
+ }
+
+
+
+ template <int dim, int spacedim>
+ typename parallel::distributed::Triangulation<dim,spacedim>::cell_iterator
+ ParticleAccessor<dim,spacedim>::get_surrounding_cell (const parallel::distributed::Triangulation<dim,spacedim> &triangulation) const
+ {
+ Assert(particle != map->end(),
+ ExcInternalError());
+
+ const typename parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell (&triangulation,
+ particle->first.first,
+ particle->first.second);
+ return cell;
+ }
+
+
+
+ template <int dim, int spacedim>
+ const ArrayView<double>
+ ParticleAccessor<dim,spacedim>::get_properties ()
+ {
+ Assert(particle != map->end(),
+ ExcInternalError());
+
+ return particle->second.get_properties();
+ }
+
+
+
+ template <int dim, int spacedim>
+ std::size_t
+ ParticleAccessor<dim,spacedim>::size () const
+ {
+ Assert(particle != map->end(),
+ ExcInternalError());
+
+ return particle->second.size();
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ ParticleAccessor<dim,spacedim>::next ()
+ {
+ Assert (particle != map->end(),ExcInternalError());
+ ++particle;
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ ParticleAccessor<dim,spacedim>::prev ()
+ {
+ Assert (particle != map->begin(),ExcInternalError());
+ --particle;
+ }
+
+
+
+ template <int dim, int spacedim>
+ bool
+ ParticleAccessor<dim,spacedim>::operator != (const ParticleAccessor<dim,spacedim> &other) const
+ {
+ return (map != other.map) || (particle != other.particle);
+ }
+
+
+
+ template <int dim, int spacedim>
+ bool
+ ParticleAccessor<dim,spacedim>::operator == (const ParticleAccessor<dim,spacedim> &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
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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 <deal_II_dimension,deal_II_space_dimension>;
+ \}
+#endif
+}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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/particles/particle_iterator.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace Particles
+{
+ template <int dim, int spacedim>
+ ParticleIterator<dim,spacedim>::ParticleIterator (const std::multimap<types::LevelInd, Particle<dim,spacedim> > &map,
+ const typename std::multimap<types::LevelInd, Particle<dim,spacedim> >::iterator &particle)
+ :
+ accessor (map, particle)
+ {}
+
+
+
+ template <int dim, int spacedim>
+ ParticleAccessor<dim,spacedim> &
+ ParticleIterator<dim,spacedim>::operator *()
+ {
+ return accessor;
+ }
+
+
+
+ template <int dim, int spacedim>
+ ParticleAccessor<dim,spacedim> *
+ ParticleIterator<dim,spacedim>::operator ->()
+ {
+ return &(this->operator* ());
+ }
+
+
+
+ template <int dim, int spacedim>
+ const ParticleAccessor<dim,spacedim> &
+ ParticleIterator<dim,spacedim>::operator *() const
+ {
+ return accessor;
+ }
+
+
+
+ template <int dim, int spacedim>
+ const ParticleAccessor<dim,spacedim> *
+ ParticleIterator<dim,spacedim>::operator ->() const
+ {
+ return &(this->operator* ());
+ }
+
+
+
+ template <int dim, int spacedim>
+ bool
+ ParticleIterator<dim,spacedim>::operator != (const ParticleIterator<dim,spacedim> &other) const
+ {
+ return accessor != other.accessor;
+ }
+
+
+
+ template <int dim, int spacedim>
+ bool
+ ParticleIterator<dim,spacedim>::operator == (const ParticleIterator<dim,spacedim> &other) const
+ {
+ return accessor == other.accessor;
+ }
+
+
+
+ template <int dim, int spacedim>
+ ParticleIterator<dim,spacedim> &
+ ParticleIterator<dim,spacedim>::operator++()
+ {
+ accessor.next();
+ return *this;
+ }
+
+
+
+ template <int dim, int spacedim>
+ ParticleIterator<dim,spacedim>
+ ParticleIterator<dim,spacedim>::operator++(int)
+ {
+ ParticleIterator tmp(*this);
+ operator++ ();
+
+ return tmp;
+ }
+
+
+
+ template <int dim, int spacedim>
+ ParticleIterator<dim,spacedim> &
+ ParticleIterator<dim,spacedim>::operator--()
+ {
+ accessor.prev();
+ return *this;
+ }
+
+
+
+ template <int dim, int spacedim>
+ ParticleIterator<dim,spacedim>
+ ParticleIterator<dim,spacedim>::operator--(int)
+ {
+ ParticleIterator tmp(*this);
+ operator-- ();
+
+ return tmp;
+ }
+}
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+DEAL_II_NAMESPACE_OPEN
+
+#include "particle_iterator.inst"
+
+DEAL_II_NAMESPACE_CLOSE
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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 <deal_II_dimension,deal_II_space_dimension>;
+ \}
+#endif
+}