--- /dev/null
+New: There is now support for the storage of Particles and their
+properties in the new namespace Particles.
+<br>
+(Rene Gassmoeller, 2017/08/09)
#include <deal.II/base/types.h>
#include <deal.II/base/array_view.h>
-#include <boost/serialization/vector.hpp>
-
DEAL_II_NAMESPACE_OPEN
/**
* A namespace that contains all classes that are related to the particle
- * implementation, in particular the fundamental <code>Particle</code> class.
+ * implementation, in particular the fundamental Particle class.
*/
namespace Particles
{
/**
* Base class of particles - represents a particle with position,
* an ID number and a variable number of properties. This class
- * can be extended to include data related to a particle by the property
- * manager.
+ * can be extended to include arbitrary data related to a particle by
+ * attaching a property pool class to the particle.
*
*
* @ingroup Particle
+ * @author Rene Gassmoeller, 2017
*
*/
- template <int dim>
+ template <int dim, int spacedim=dim, typename PropertyType = double>
class Particle
{
public:
/**
* Constructor for Particle, creates a particle with the specified
- * ID at the specified location. Note that deal.II
- * does not check for duplicate particle IDs so the user must
+ * ID at the specified location. Note that there is no
+ * check for duplicate particle IDs so the user must
* make sure the IDs are unique over all processes.
*
- * @param new_location Initial location of particle.
- * @param new_reference_location Initial location of the particle
+ * @param[in] location Initial location of particle.
+ * @param[in] reference_location Initial location of the particle
* in the coordinate system of the reference cell.
- * @param new_id Globally unique ID number of particle.
+ * @param[in] id Globally unique ID number of particle.
*/
- Particle (const Point<dim> &new_location,
- const Point<dim> &new_reference_location,
- const types::particle_index new_id);
+ Particle (const Point<spacedim> &location,
+ const Point<dim> &reference_location,
+ const types::particle_index id);
/**
* Copy-Constructor for Particle, creates a particle with exactly the
* for registering and freeing this memory in the property pool this
* constructor registers a new chunk, and copies the properties.
*/
- Particle (const Particle<dim> &particle);
+ Particle (const Particle<dim,spacedim,PropertyType> &particle);
/**
* Constructor for Particle, creates a particle from a data vector.
- * This constructor is usually called after sending a particle to a
- * different process.
+ * This constructor is usually called after serializing a particle by
+ * calling the write_data function.
*
* @param[in,out] begin_data A pointer to a memory location from which
* to read the information that completely describes a particle. This
* allocate the property data used by this particle.
*/
Particle (const void *&begin_data,
- PropertyPool &new_property_pool);
+ PropertyPool<PropertyType> &new_property_pool);
/**
* Move constructor for Particle, creates a particle from an existing
- * one by stealing its state.
+ * one by stealing its state and properties.
*/
- Particle (Particle<dim> &&particle);
+ Particle (Particle<dim,spacedim,PropertyType> &&particle);
/**
* Copy assignment operator.
*/
- Particle<dim> &operator=(const Particle<dim> &particle);
+ Particle<dim,spacedim,PropertyType> &operator=(const Particle<dim,spacedim,PropertyType> &particle);
/**
* Move assignment operator.
*/
- Particle<dim> &operator=(Particle<dim> &&particle);
+ Particle<dim,spacedim,PropertyType> &operator=(Particle<dim,spacedim,PropertyType> &&particle);
/**
* Destructor. Releases the property handle if it is valid, and
* therefore frees that memory space for other particles. (Note:
- * the memory is managed by the property_pool, and the memory is not
- * deallocated by this function, it is kept in reserve for other
- * particles).
+ * the memory is managed by the property pool, and the pool is responsible
+ * for what happens to the memory.
*/
~Particle ();
* 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_loc The new location for this particle.
+ * @param [in] new_location The new location for this particle.
*/
void
- set_location (const Point<dim> &new_loc);
+ set_location (const Point<spacedim> &new_location);
/**
* Get the location of this particle.
*
* @return The location of this particle.
*/
- const Point<dim> &
+ const Point<spacedim> &
get_location () const;
/**
* Set the reference location of this particle.
*
- * @param [in] new_loc The new reference location for this particle.
+ * @param [in] new_reference_location The new reference location for
+ * this particle.
*/
void
- set_reference_location (const Point<dim> &new_loc);
+ set_reference_location (const Point<dim> &new_reference_location);
/**
* Return the reference location of this particle in its current cell.
/**
* 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
+ * since the particle does not know about the properties
+ * we do not want to do it at construction time. Another use for this
* function is after particle transfer to a new process.
*/
void
- set_property_pool(PropertyPool &property_pool);
+ set_property_pool(PropertyPool<PropertyType> &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.
+ * @param [in] new_properties An ArrayView of the new
+ * properties that will be copied into the property pool.
*/
void
- set_properties (const std::vector<double> &new_properties);
+ set_properties (const ArrayView<const PropertyType> &new_properties);
/**
* Get write-access to properties of this particle.
*
* @return An ArrayView of the properties of this particle.
*/
- const ArrayView<double>
+ const ArrayView<PropertyType>
get_properties ();
/**
*
* @return An ArrayView of the properties of this particle.
*/
- const ArrayView<const double>
+ const ArrayView<const PropertyType>
get_properties () const;
/**
/**
* Current particle location.
*/
- Point<dim> location;
+ Point<spacedim> location;
/**
* Current particle location in the reference cell.
* A pointer to the property pool. Necessary to translate from the
* handle to the actual memory locations.
*/
- PropertyPool *property_pool;
+ PropertyPool<PropertyType> *property_pool;
/**
- * A handle to all particle properties
+ * A handle to all particle properties.
*/
- PropertyPool::Handle properties;
+ typename PropertyPool<PropertyType>::Handle properties;
};
/* -------------------------- inline and template functions ---------------------- */
- template <int dim>
+ template <int dim, int spacedim, typename PropertyType>
template <class Archive>
- void Particle<dim>::serialize (Archive &ar, const unsigned int)
+ void Particle<dim,spacedim,PropertyType>::serialize (Archive &ar, const unsigned int)
{
ar &location
& id
namespace Particles
{
/**
- * This class manages the memory space in which particles store their
- * properties. Because this is dynamic memory and every particle needs the
- * same amount it is more efficient to let this be handled by a central
- * manager that does not need to allocate/deallocate memory every time a
- * particle is constructed/destroyed.
+ * This class manages a memory space in which particles store their
+ * properties. Because this is dynamic memory and often every particle
+ * needs the same amount, it is more efficient to let this be handled by a
+ * central manager that does not need to allocate/deallocate memory every
+ * time a particle is constructed/destroyed.
+ * The current implementation uses simple new/delete allocation for every
+ * block. Additionally, the current implementation
+ * assumes the same number of properties per particle, but of course the
+ * PropertyType could contain a pointer to dynamically allocated memory
+ * with varying sizes per particle (this memory would not be managed by this
+ * class).
+ * Because PropertyPool only returns handles it could be enhanced internally
+ * (e.g. to allow for varying number of properties per handle) without
+ * affecting its interface.
*/
+ template <typename PropertyType = double>
class PropertyPool
{
public:
* uniquely identifies the slot of memory that is reserved for this
* particle.
*/
- typedef double *Handle;
+ typedef PropertyType *Handle;
/**
* Define a default (invalid) value for handles.
/**
* Constructor. Stores the number of properties per reserved slot.
+ * By default this class assumes the template parameter PropertyType
+ * completely describes all properties of the particle.
*/
- PropertyPool (const unsigned int n_properties_per_slot);
+ PropertyPool (const unsigned int n_properties_per_slot = 1);
/**
* Returns a new handle that allows accessing the reserved block
- * of memory.
+ * of memory. The returned memory block is not initialized.
*/
Handle allocate_properties_array ();
/**
* Mark the properties corresponding to the handle @p handle as
- * deleted. After calling this function calling get_properties() with
- * the freed handle causes undefined behavior.
+ * deleted. After calling this function @p handle is invalid and
+ * can not longer be used to access properties.
*/
- void deallocate_properties_array (const Handle handle);
+ void deallocate_properties_array (Handle &handle);
/**
* Return an ArrayView to the properties that correspond to the given
- * handle @p handle.
+ * handle @p handle. This function allows read-write access to the
+ * underlying properties.
*/
- ArrayView<double> get_properties (const Handle handle);
+ const ArrayView<PropertyType> get_properties (const Handle handle);
/**
- * Reserves the dynamic memory needed for storing the properties of
- * @p size particles.
+ * Return an ArrayView to the properties that correspond to the given
+ * handle @p handle. This function only allows read access to the underlying
+ * properties.
*/
- void reserve(const std::size_t size);
+ const ArrayView<const PropertyType> get_properties (const Handle handle) const;
private:
/**
--- /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_property_pool_templates_h
+#define dealii__particles_property_pool_templates_h
+
+#include <deal.II/particles/property_pool.h>
+#include <deal.II/base/exceptions.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace Particles
+{
+ template <typename PropertyType>
+ const typename PropertyPool<PropertyType>::Handle PropertyPool<PropertyType>::invalid_handle = NULL;
+
+
+
+ template <typename PropertyType>
+ PropertyPool<PropertyType>::PropertyPool (const unsigned int n_properties_per_slot)
+ :
+ n_properties (n_properties_per_slot)
+ {}
+
+
+
+ template <typename PropertyType>
+ typename PropertyPool<PropertyType>::Handle
+ PropertyPool<PropertyType>::allocate_properties_array ()
+ {
+ return new PropertyType[n_properties];
+ }
+
+
+ template <typename PropertyType>
+ void
+ PropertyPool<PropertyType>::deallocate_properties_array (Handle &handle)
+ {
+ delete[] handle;
+ handle = invalid_handle;
+ }
+
+
+
+ template <typename PropertyType>
+ const ArrayView<PropertyType>
+ PropertyPool<PropertyType>::get_properties (const Handle handle)
+ {
+ Assert(handle != PropertyPool<PropertyType>::invalid_handle,
+ ExcMessage("PropertyPool was asked for the properties that belong to an invalid handle. "
+ "Either this handle was never validly created, or it was deallocated before use."));
+
+ return ArrayView<PropertyType>(handle, n_properties);
+ }
+
+
+
+ template <typename PropertyType>
+ const ArrayView<const PropertyType>
+ PropertyPool<PropertyType>::get_properties (const Handle handle) const
+ {
+ Assert(handle != PropertyPool<PropertyType>::invalid_handle,
+ ExcMessage("PropertyPool was asked for the properties that belong to an invalid handle. "
+ "Either this handle was never validly created, or it was deallocated before use."));
+
+ return ArrayView<const PropertyType>(handle, n_properties);
+ }
+}
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
SET(_inst
particle.inst.in
+ property_pool.inst.in
)
FILE(GLOB _header
namespace Particles
{
- template <int dim>
- Particle<dim>::Particle ()
+ template <int dim, int spacedim, typename PropertyType>
+ Particle<dim,spacedim,PropertyType>::Particle ()
:
location (),
reference_location(),
id (0),
property_pool(NULL),
- properties(PropertyPool::invalid_handle)
+ properties(PropertyPool<PropertyType>::invalid_handle)
{
}
- template <int dim>
- Particle<dim>::Particle (const Point<dim> &new_location,
- const Point<dim> &new_reference_location,
- const types::particle_index new_id)
+
+ template <int dim, int spacedim, typename PropertyType>
+ Particle<dim,spacedim,PropertyType>::Particle (const Point<spacedim> &location,
+ const Point<dim> &reference_location,
+ const types::particle_index id)
:
- location (new_location),
- reference_location (new_reference_location),
- id (new_id),
+ location (location),
+ reference_location (reference_location),
+ id (id),
property_pool(NULL),
- properties (PropertyPool::invalid_handle)
+ properties (PropertyPool<PropertyType>::invalid_handle)
{
}
- template <int dim>
- Particle<dim>::Particle (const Particle<dim> &particle)
+
+
+ template <int dim, int spacedim, typename PropertyType>
+ Particle<dim,spacedim,PropertyType>::Particle (const Particle<dim,spacedim,PropertyType> &particle)
:
location (particle.get_location()),
reference_location (particle.get_reference_location()),
id (particle.get_id()),
property_pool(particle.property_pool),
- properties ((property_pool != NULL) ? property_pool->allocate_properties_array() : PropertyPool::invalid_handle)
+ properties ((property_pool != NULL) ? property_pool->allocate_properties_array() : PropertyPool<PropertyType>::invalid_handle)
{
if (property_pool != NULL)
{
- const ArrayView<double> my_properties = property_pool->get_properties(properties);
+ const ArrayView<PropertyType> my_properties = property_pool->get_properties(properties);
if (my_properties.size() != 0)
{
- const ArrayView<const double> their_properties = particle.get_properties();
+ const ArrayView<const PropertyType> their_properties = particle.get_properties();
std::copy(&their_properties[0],&their_properties[0]+their_properties.size(),&my_properties[0]);
}
}
- template <int dim>
- Particle<dim>::Particle (const void *&data,
- PropertyPool &new_property_pool)
+
+ template <int dim, int spacedim, typename PropertyType>
+ Particle<dim,spacedim,PropertyType>::Particle (const void *&data,
+ PropertyPool<PropertyType> &new_property_pool)
{
const types::particle_index *id_data = static_cast<const types::particle_index *> (data);
id = *id_data++;
property_pool = &new_property_pool;
properties = property_pool->allocate_properties_array();
+ // TODO reinterpret cast pdata to PropertyType
// See if there are properties to load
- const ArrayView<double> particle_properties = property_pool->get_properties(properties);
+ const ArrayView<PropertyType> particle_properties = property_pool->get_properties(properties);
for (unsigned int i = 0; i < particle_properties.size(); ++i)
particle_properties[i] = *pdata++;
data = static_cast<const void *> (pdata);
}
- template <int dim>
- Particle<dim>::Particle (Particle<dim> &&particle)
+
+
+ template <int dim, int spacedim, typename PropertyType>
+ Particle<dim,spacedim,PropertyType>::Particle (Particle<dim,spacedim,PropertyType> &&particle)
:
location (particle.location),
reference_location(particle.reference_location),
property_pool(particle.property_pool),
properties(particle.properties)
{
- particle.properties = PropertyPool::invalid_handle;
+ particle.properties = PropertyPool<PropertyType>::invalid_handle;
}
- template <int dim>
- Particle<dim> &
- Particle<dim>::operator=(const Particle<dim> &particle)
+
+
+ template <int dim, int spacedim, typename PropertyType>
+ Particle<dim,spacedim,PropertyType> &
+ Particle<dim,spacedim,PropertyType>::operator=(const Particle<dim,spacedim,PropertyType> &particle)
{
if (this != &particle)
{
if (property_pool != NULL)
{
properties = property_pool->allocate_properties_array();
- const ArrayView<const double> their_properties = particle.get_properties();
+ const ArrayView<const PropertyType> their_properties = particle.get_properties();
if (their_properties.size() != 0)
{
- const ArrayView<double> my_properties = property_pool->get_properties(properties);
+ const ArrayView<PropertyType> my_properties = property_pool->get_properties(properties);
std::copy(&their_properties[0],&their_properties[0]+their_properties.size(),&my_properties[0]);
}
}
else
- properties = PropertyPool::invalid_handle;
+ properties = PropertyPool<PropertyType>::invalid_handle;
}
return *this;
}
- template <int dim>
- Particle<dim> &
- Particle<dim>::operator=(Particle<dim> &&particle)
+
+
+ template <int dim, int spacedim, typename PropertyType>
+ Particle<dim,spacedim,PropertyType> &
+ Particle<dim,spacedim,PropertyType>::operator=(Particle<dim,spacedim,PropertyType> &&particle)
{
if (this != &particle)
{
id = particle.id;
property_pool = particle.property_pool;
properties = particle.properties;
- particle.properties = PropertyPool::invalid_handle;
+ particle.properties = PropertyPool<PropertyType>::invalid_handle;
}
return *this;
}
- template <int dim>
- Particle<dim>::~Particle ()
+
+
+ template <int dim, int spacedim, typename PropertyType>
+ Particle<dim,spacedim,PropertyType>::~Particle ()
{
- if (properties != PropertyPool::invalid_handle)
+ if (properties != PropertyPool<PropertyType>::invalid_handle)
property_pool->deallocate_properties_array(properties);
}
- template <int dim>
+
+
+ template <int dim, int spacedim, typename PropertyType>
void
- Particle<dim>::write_data (void *&data) const
+ Particle<dim,spacedim,PropertyType>::write_data (void *&data) const
{
types::particle_index *id_data = static_cast<types::particle_index *> (data);
*id_data = id;
for (unsigned int i = 0; i < dim; ++i,++pdata)
*pdata = reference_location(i);
+ // TODO reinterpret cast pdata to PropertyType
// Write property data
- const ArrayView<double> particle_properties = property_pool->get_properties(properties);
+ const ArrayView<PropertyType> particle_properties = property_pool->get_properties(properties);
for (unsigned int i = 0; i < particle_properties.size(); ++i,++pdata)
*pdata = particle_properties[i];
data = static_cast<void *> (pdata);
}
- template <int dim>
+
+
+ template <int dim, int spacedim, typename PropertyType>
void
- Particle<dim>::set_location (const Point<dim> &new_loc)
+ Particle<dim,spacedim,PropertyType>::set_location (const Point<spacedim> &new_loc)
{
location = new_loc;
}
- template <int dim>
- const Point<dim> &
- Particle<dim>::get_location () const
+
+
+ template <int dim, int spacedim, typename PropertyType>
+ const Point<spacedim> &
+ Particle<dim,spacedim,PropertyType>::get_location () const
{
return location;
}
- template <int dim>
+
+
+ template <int dim, int spacedim, typename PropertyType>
void
- Particle<dim>::set_reference_location (const Point<dim> &new_loc)
+ Particle<dim,spacedim,PropertyType>::set_reference_location (const Point<dim> &new_loc)
{
reference_location = new_loc;
}
- template <int dim>
+
+
+ template <int dim, int spacedim, typename PropertyType>
const Point<dim> &
- Particle<dim>::get_reference_location () const
+ Particle<dim,spacedim,PropertyType>::get_reference_location () const
{
return reference_location;
}
- template <int dim>
+
+
+ template <int dim, int spacedim, typename PropertyType>
types::particle_index
- Particle<dim>::get_id () const
+ Particle<dim,spacedim,PropertyType>::get_id () const
{
return id;
}
- template <int dim>
+
+
+ template <int dim, int spacedim, typename PropertyType>
void
- Particle<dim>::set_property_pool (PropertyPool &new_property_pool)
+ Particle<dim,spacedim,PropertyType>::set_property_pool (PropertyPool<PropertyType> &new_property_pool)
{
property_pool = &new_property_pool;
}
- template <int dim>
+
+
+ template <int dim, int spacedim, typename PropertyType>
+ bool
+ Particle<dim,spacedim,PropertyType>::has_properties () const
+ {
+ return (property_pool != NULL)
+ && (properties != PropertyPool<PropertyType>::invalid_handle);
+ }
+
+
+
+ template <int dim, int spacedim, typename PropertyType>
void
- Particle<dim>::set_properties (const std::vector<double> &new_properties)
+ Particle<dim,spacedim,PropertyType>::set_properties (const ArrayView<const PropertyType> &new_properties)
{
- if (properties == PropertyPool::invalid_handle)
+ if (properties == PropertyPool<PropertyType>::invalid_handle)
properties = property_pool->allocate_properties_array();
- const ArrayView<double> old_properties = property_pool->get_properties(properties);
+ const ArrayView<PropertyType> old_properties = property_pool->get_properties(properties);
std::copy(new_properties.begin(),new_properties.end(),&old_properties[0]);
}
- template <int dim>
- const ArrayView<const double>
- Particle<dim>::get_properties () const
+
+
+ template <int dim, int spacedim, typename PropertyType>
+ const ArrayView<const PropertyType>
+ Particle<dim,spacedim,PropertyType>::get_properties () const
{
Assert(property_pool != NULL,
ExcInternalError());
}
- template <int dim>
- const ArrayView<double>
- Particle<dim>::get_properties ()
+
+ template <int dim, int spacedim, typename PropertyType>
+ const ArrayView<PropertyType>
+ Particle<dim,spacedim,PropertyType>::get_properties ()
{
Assert(property_pool != NULL,
ExcInternalError());
// ---------------------------------------------------------------------
-for (deal_II_dimension : DIMENSIONS)
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS; number : REAL_SCALARS)
{
namespace Particles
\{
template
- class Particle <deal_II_dimension>;
+ class Particle <deal_II_dimension,deal_II_space_dimension,number>;
\}
}
// ---------------------------------------------------------------------
-#include <deal.II/particles/property_pool.h>
+#include <deal.II/particles/property_pool.templates.h>
DEAL_II_NAMESPACE_OPEN
-namespace Particles
-{
- const typename PropertyPool::Handle PropertyPool::invalid_handle = NULL;
-
-
- PropertyPool::PropertyPool (const unsigned int n_properties_per_slot)
- :
- n_properties (n_properties_per_slot)
- {}
-
-
-
- PropertyPool::Handle
- PropertyPool::allocate_properties_array ()
- {
- return new double[n_properties];
- }
-
-
-
- void
- PropertyPool::deallocate_properties_array (Handle handle)
- {
- delete[] handle;
- }
-
-
-
- ArrayView<double>
- PropertyPool::get_properties (const Handle handle)
- {
- return ArrayView<double>(handle, n_properties);
- }
-
-
-
- void
- PropertyPool::reserve(const std::size_t size)
- {
- (void)size;
- }
-}
+#include "property_pool.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 (number : REAL_SCALARS)
+{
+ namespace Particles
+ \{
+ template
+ class PropertyPool <number>;
+ \}
+}
#include "../tests.h"
#include <deal.II/particles/particle.h>
-#include <fstream>
-#include <iomanip>
template <int dim>
#include "../tests.h"
#include <deal.II/particles/particle.h>
-#include <fstream>
-#include <iomanip>
template <int dim>
Particles::Particle<2> particle(position,reference_position,index);
deallog << "Particle location: " << particle.get_location() << std::endl
- << "Particle reference location: " << particle.get_reference_location() << std::endl
- << "Particle index: " << particle.get_id() << std::endl;
+ << "Particle reference location: " << particle.get_reference_location() << std::endl
+ << "Particle index: " << particle.get_id() << std::endl;
const Particles::Particle<2> copy(particle);
deallog << "Copy particle location: " << copy.get_location() << std::endl
- << "Copy particle reference location: " << copy.get_reference_location() << std::endl
- << "Copy particle index: " << copy.get_id() << std::endl;
+ << "Copy particle reference location: " << copy.get_reference_location() << std::endl
+ << "Copy particle index: " << copy.get_id() << std::endl;
const Particles::Particle<2> moved_particle(std::move(particle));
deallog << "Moved particle location: " << moved_particle.get_location() << std::endl
- << "Moved particle reference location: " << moved_particle.get_reference_location() << std::endl
- << "Moved particle index: " << moved_particle.get_id() << std::endl;
+ << "Moved particle reference location: " << moved_particle.get_reference_location() << std::endl
+ << "Moved particle index: " << moved_particle.get_id() << std::endl;
deallog << "Original particle location: " << particle.get_location() << std::endl
- << "Original particle reference location: " << particle.get_reference_location() << std::endl
- << "Original particle index: " << particle.get_id() << std::endl;
+ << "Original particle reference location: " << particle.get_reference_location() << std::endl
+ << "Original particle index: " << particle.get_id() << std::endl;
}
deallog << "OK" << std::endl;
#include "../tests.h"
#include <deal.II/particles/particle.h>
-#include <fstream>
-#include <iomanip>
+#include <deal.II/base/array_view.h>
template <int dim>
{
{
const unsigned int n_properties_per_particle = 3;
- Particles::PropertyPool pool(n_properties_per_particle);
+ Particles::PropertyPool<> pool(n_properties_per_particle);
Point<2> position;
position(0) = 0.3;
Particles::Particle<2> particle(position,reference_position,index);
particle.set_property_pool(pool);
- particle.set_properties(properties);
+ particle.set_properties(ArrayView<double>(&properties[0],properties.size()));
deallog << "Particle properties: "
- << std::vector<double>(particle.get_properties().begin(),particle.get_properties().end())
- << std::endl;
+ << std::vector<double>(particle.get_properties().begin(),particle.get_properties().end())
+ << std::endl;
const Particles::Particle<2> copy(particle);
deallog << "Copy particle properties: "
- << std::vector<double>(copy.get_properties().begin(),copy.get_properties().end())
- << std::endl;
+ << std::vector<double>(copy.get_properties().begin(),copy.get_properties().end())
+ << std::endl;
+
+ deallog << "Old particle has properties before move: " << particle.has_properties() << std::endl;
const Particles::Particle<2> moved_particle(std::move(particle));
+ deallog << "Old particle has properties after move: " << particle.has_properties() << std::endl;
+
deallog << "Moved particle properties: "
- << std::vector<double>(moved_particle.get_properties().begin(),moved_particle.get_properties().end())
- << std::endl;
+ << std::vector<double>(moved_particle.get_properties().begin(),moved_particle.get_properties().end())
+ << std::endl;
+
}
deallog << "OK" << std::endl;
DEAL::Particle properties: 0.150000 0.450000 0.750000
DEAL::Copy particle properties: 0.150000 0.450000 0.750000
+DEAL::Old particle has properties before move: 1
+DEAL::Old particle has properties after move: 0
DEAL::Moved particle properties: 0.150000 0.450000 0.750000
DEAL::OK
--- /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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// check the creation, simplest usage, and destruction of a property pool
+
+#include "../tests.h"
+#include <deal.II/particles/property_pool.h>
+#include <fstream>
+#include <iomanip>
+
+
+template <typename PropertyType>
+void test ()
+{
+ {
+ Particles::PropertyPool<PropertyType> pool;
+
+ typename Particles::PropertyPool<PropertyType>::Handle handle = pool.allocate_properties_array();
+
+ pool.get_properties(handle)[0] = 2.5;
+
+ deallog << "Pool properties: " << pool.get_properties(handle)[0] << std::endl;
+
+ pool.deallocate_properties_array(handle);
+ }
+
+ deallog << "OK" << std::endl;
+}
+
+
+
+int main ()
+{
+ initlog();
+ test<double>();
+ test<float>();
+}
--- /dev/null
+
+DEAL::Pool properties: 2.50000
+DEAL::OK
+DEAL::Pool properties: 2.50000
+DEAL::OK
--- /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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// test a property pool that allocates more than one property per chunk
+
+#include "../tests.h"
+#include <deal.II/particles/property_pool.h>
+#include <fstream>
+#include <iomanip>
+
+
+template <typename PropertyType>
+void test ()
+{
+ {
+ const unsigned int n_properties = 3;
+ Particles::PropertyPool<PropertyType> pool(n_properties);
+
+ typename Particles::PropertyPool<PropertyType>::Handle handle = pool.allocate_properties_array();
+
+ pool.get_properties(handle)[0] = 1.2;
+ pool.get_properties(handle)[1] = 2.5;
+ pool.get_properties(handle)[2] = 2.7;
+
+
+ deallog << "Pool properties:";
+
+ for (unsigned int i=0; i<pool.get_properties(handle).size(); ++i)
+ deallog << " " << pool.get_properties(handle)[i];
+
+ deallog << std::endl;
+
+ pool.deallocate_properties_array(handle);
+ }
+
+ deallog << "OK" << std::endl;
+}
+
+
+
+int main ()
+{
+ initlog();
+ test<double>();
+ test<float>();
+}
--- /dev/null
+
+DEAL::Pool properties: 1.20000 2.50000 2.70000
+DEAL::OK
+DEAL::Pool properties: 1.20000 2.50000 2.70000
+DEAL::OK
--- /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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// test a property pool that uses a custom type as storage container
+
+#include "../tests.h"
+#include <deal.II/particles/property_pool.templates.h>
+#include <fstream>
+#include <iomanip>
+
+
+
+struct DataType
+{
+ double a[2];
+ char c;
+ double *pointer;
+};
+
+void test ()
+{
+ {
+ double b[2] = {0.7,1.3};
+
+ Particles::PropertyPool<DataType> pool;
+
+ typename Particles::PropertyPool<DataType>::Handle handle = pool.allocate_properties_array();
+
+ DataType &properties_in = pool.get_properties(handle)[0];
+ properties_in.a[0] = 2.5;
+ properties_in.a[1] = 2.0;
+ properties_in.c = 'f';
+ properties_in.pointer = &b[1];
+
+ const DataType &properties_out = pool.get_properties(handle)[0];
+
+ deallog << "Pool properties: " << properties_out.a[0] << ' '
+ << properties_out.a[1] << ' '
+ << properties_out.c << ' '
+ << *properties_out.pointer
+ << std::endl;
+
+ pool.deallocate_properties_array(handle);
+ }
+
+ deallog << "OK" << std::endl;
+}
+
+
+
+int main ()
+{
+ initlog();
+ ;
+
+ test();
+}
--- /dev/null
+
+DEAL::Pool properties: 2.50000 2.00000 f 1.30000
+DEAL::OK