]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add accessor and iterator classes
authorRene Gassmoeller <rene.gassmoeller@mailbox.org>
Wed, 13 Sep 2017 15:18:10 +0000 (17:18 +0200)
committerRene Gassmoeller <rene.gassmoeller@mailbox.org>
Fri, 22 Sep 2017 15:59:07 +0000 (09:59 -0600)
include/deal.II/particles/particle.h
include/deal.II/particles/particle_accessor.h [new file with mode: 0644]
include/deal.II/particles/particle_iterator.h [new file with mode: 0644]
source/particles/CMakeLists.txt
source/particles/particle.cc
source/particles/particle.inst.in
source/particles/particle_accessor.cc [new file with mode: 0644]
source/particles/particle_accessor.inst.in [new file with mode: 0644]
source/particles/particle_iterator.cc [new file with mode: 0644]
source/particles/particle_iterator.inst.in [new file with mode: 0644]

index 9a2de57d1120ea896ddbbc60f668c2cebcc11714..d27082cb92cc205a39099ae594f9212dc5c00063 100644 (file)
@@ -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 <deal.II/particles/property_pool.h>
 
@@ -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<int, int> 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 (file)
index 0000000..94c5142
--- /dev/null
@@ -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 <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
diff --git a/include/deal.II/particles/particle_iterator.h b/include/deal.II/particles/particle_iterator.h
new file mode 100644 (file)
index 0000000..83941dd
--- /dev/null
@@ -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/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
index f355e5050f4e906b182b14ce78c664065a9502be..8eaa91767d1eb5a1535e309edcc9f52d5b70e49b 100644 (file)
@@ -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
   )
 
index e2917c5f00d65f40ab509d50ea52f54c2bc04dc9..b898f21193f8930d26257b75633d8f6de20b1e52 100644 (file)
@@ -215,8 +215,8 @@ namespace Particles
   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())
       {
@@ -287,7 +287,7 @@ namespace Particles
   Particle<dim,spacedim,PropertyType>::has_properties () const
   {
     return (property_pool != NULL)
-        && (properties != PropertyPool<PropertyType>::invalid_handle);
+           && (properties != PropertyPool<PropertyType>::invalid_handle);
   }
 
 
@@ -297,8 +297,8 @@ namespace Particles
   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();
@@ -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);
   }
index 9ea1b7927dc999552842830bfaf968c2507f6e57..1f1327202ada708e3577ea79b3b756c0597d648a 100644 (file)
 
 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
 }
diff --git a/source/particles/particle_accessor.cc b/source/particles/particle_accessor.cc
new file mode 100644 (file)
index 0000000..16bbd14
--- /dev/null
@@ -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/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
diff --git a/source/particles/particle_accessor.inst.in b/source/particles/particle_accessor.inst.in
new file mode 100644 (file)
index 0000000..71a9286
--- /dev/null
@@ -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 <deal_II_dimension,deal_II_space_dimension>;
+    \}
+#endif
+}
diff --git a/source/particles/particle_iterator.cc b/source/particles/particle_iterator.cc
new file mode 100644 (file)
index 0000000..4fa7873
--- /dev/null
@@ -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/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
diff --git a/source/particles/particle_iterator.inst.in b/source/particles/particle_iterator.inst.in
new file mode 100644 (file)
index 0000000..2d6e7c9
--- /dev/null
@@ -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 <deal_II_dimension,deal_II_space_dimension>;
+    \}
+#endif
+}

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.