From f245d923158520866971c534ea2d15749b923566 Mon Sep 17 00:00:00 2001 From: Rene Gassmoeller Date: Wed, 9 Aug 2017 12:12:07 -0600 Subject: [PATCH] Templatize property pool and particle. Address comments --- .../changes/major/20170809ReneGassmoeller | 4 + include/deal.II/particles/particle.h | 94 +++++----- include/deal.II/particles/property_pool.h | 46 +++-- .../particles/property_pool.templates.h | 84 +++++++++ source/particles/CMakeLists.txt | 1 + source/particles/particle.cc | 167 +++++++++++------- source/particles/particle.inst.in | 4 +- source/particles/property_pool.cc | 45 +---- source/particles/property_pool.inst.in | 24 +++ tests/particles/particle_01.cc | 2 - tests/particles/particle_02.cc | 18 +- tests/particles/particle_03.cc | 24 +-- tests/particles/particle_03.output | 2 + tests/particles/property_pool_01.cc | 51 ++++++ tests/particles/property_pool_01.output | 5 + tests/particles/property_pool_02.cc | 60 +++++++ tests/particles/property_pool_02.output | 5 + tests/particles/property_pool_03.cc | 71 ++++++++ tests/particles/property_pool_03.output | 3 + 19 files changed, 521 insertions(+), 189 deletions(-) create mode 100644 doc/news/changes/major/20170809ReneGassmoeller create mode 100644 include/deal.II/particles/property_pool.templates.h create mode 100644 source/particles/property_pool.inst.in create mode 100644 tests/particles/property_pool_01.cc create mode 100644 tests/particles/property_pool_01.output create mode 100644 tests/particles/property_pool_02.cc create mode 100644 tests/particles/property_pool_02.output create mode 100644 tests/particles/property_pool_03.cc create mode 100644 tests/particles/property_pool_03.output diff --git a/doc/news/changes/major/20170809ReneGassmoeller b/doc/news/changes/major/20170809ReneGassmoeller new file mode 100644 index 0000000000..354eda6d6a --- /dev/null +++ b/doc/news/changes/major/20170809ReneGassmoeller @@ -0,0 +1,4 @@ +New: There is now support for the storage of Particles and their +properties in the new namespace Particles. +
+(Rene Gassmoeller, 2017/08/09) diff --git a/include/deal.II/particles/particle.h b/include/deal.II/particles/particle.h index 6c7bc369fd..cedf7c4055 100644 --- a/include/deal.II/particles/particle.h +++ b/include/deal.II/particles/particle.h @@ -22,13 +22,11 @@ #include #include -#include - DEAL_II_NAMESPACE_OPEN /** * A namespace that contains all classes that are related to the particle - * implementation, in particular the fundamental Particle class. + * implementation, in particular the fundamental Particle class. */ namespace Particles { @@ -83,14 +81,15 @@ 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 + template class Particle { public: @@ -102,18 +101,18 @@ namespace Particles /** * 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 &new_location, - const Point &new_reference_location, - const types::particle_index new_id); + Particle (const Point &location, + const Point &reference_location, + const types::particle_index id); /** * Copy-Constructor for Particle, creates a particle with exactly the @@ -122,12 +121,12 @@ namespace Particles * for registering and freeing this memory in the property pool this * constructor registers a new chunk, and copies the properties. */ - Particle (const Particle &particle); + Particle (const Particle &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 @@ -138,30 +137,29 @@ namespace Particles * allocate the property data used by this particle. */ Particle (const void *&begin_data, - PropertyPool &new_property_pool); + PropertyPool &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 &&particle); + Particle (Particle &&particle); /** * Copy assignment operator. */ - Particle &operator=(const Particle &particle); + Particle &operator=(const Particle &particle); /** * Move assignment operator. */ - Particle &operator=(Particle &&particle); + Particle &operator=(Particle &&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 (); @@ -185,26 +183,27 @@ namespace Particles * 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 &new_loc); + set_location (const Point &new_location); /** * Get the location of this particle. * * @return The location of this particle. */ - const Point & + const Point & 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 &new_loc); + set_reference_location (const Point &new_reference_location); /** * Return the reference location of this particle in its current cell. @@ -221,28 +220,35 @@ namespace Particles /** * 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 &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 &new_properties); + set_properties (const ArrayView &new_properties); /** * Get write-access to properties of this particle. * * @return An ArrayView of the properties of this particle. */ - const ArrayView + const ArrayView get_properties (); /** @@ -250,7 +256,7 @@ namespace Particles * * @return An ArrayView of the properties of this particle. */ - const ArrayView + const ArrayView get_properties () const; /** @@ -263,7 +269,7 @@ namespace Particles /** * Current particle location. */ - Point location; + Point location; /** * Current particle location in the reference cell. @@ -279,19 +285,19 @@ namespace Particles * A pointer to the property pool. Necessary to translate from the * handle to the actual memory locations. */ - PropertyPool *property_pool; + PropertyPool *property_pool; /** - * A handle to all particle properties + * A handle to all particle properties. */ - PropertyPool::Handle properties; + typename PropertyPool::Handle properties; }; /* -------------------------- inline and template functions ---------------------- */ - template + template template - void Particle::serialize (Archive &ar, const unsigned int) + void Particle::serialize (Archive &ar, const unsigned int) { ar &location & id diff --git a/include/deal.II/particles/property_pool.h b/include/deal.II/particles/property_pool.h index a392c23449..176f1bfb0b 100644 --- a/include/deal.II/particles/property_pool.h +++ b/include/deal.II/particles/property_pool.h @@ -23,12 +23,22 @@ DEAL_II_NAMESPACE_OPEN 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 class PropertyPool { public: @@ -37,7 +47,7 @@ namespace Particles * 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. @@ -46,33 +56,37 @@ namespace Particles /** * 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 get_properties (const Handle handle); + const ArrayView 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 get_properties (const Handle handle) const; private: /** diff --git a/include/deal.II/particles/property_pool.templates.h b/include/deal.II/particles/property_pool.templates.h new file mode 100644 index 0000000000..7d76c9c0e1 --- /dev/null +++ b/include/deal.II/particles/property_pool.templates.h @@ -0,0 +1,84 @@ +// --------------------------------------------------------------------- +// +// 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 +#include + +DEAL_II_NAMESPACE_OPEN + +namespace Particles +{ + template + const typename PropertyPool::Handle PropertyPool::invalid_handle = NULL; + + + + template + PropertyPool::PropertyPool (const unsigned int n_properties_per_slot) + : + n_properties (n_properties_per_slot) + {} + + + + template + typename PropertyPool::Handle + PropertyPool::allocate_properties_array () + { + return new PropertyType[n_properties]; + } + + + template + void + PropertyPool::deallocate_properties_array (Handle &handle) + { + delete[] handle; + handle = invalid_handle; + } + + + + template + const ArrayView + PropertyPool::get_properties (const Handle handle) + { + Assert(handle != PropertyPool::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(handle, n_properties); + } + + + + template + const ArrayView + PropertyPool::get_properties (const Handle handle) const + { + Assert(handle != PropertyPool::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(handle, n_properties); + } +} + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/source/particles/CMakeLists.txt b/source/particles/CMakeLists.txt index b45a281df6..f355e5050f 100644 --- a/source/particles/CMakeLists.txt +++ b/source/particles/CMakeLists.txt @@ -23,6 +23,7 @@ SET(_src SET(_inst particle.inst.in + property_pool.inst.in ) FILE(GLOB _header diff --git a/source/particles/particle.cc b/source/particles/particle.cc index 11e74a23ca..9ea08746a3 100644 --- a/source/particles/particle.cc +++ b/source/particles/particle.cc @@ -19,47 +19,50 @@ DEAL_II_NAMESPACE_OPEN namespace Particles { - template - Particle::Particle () + template + Particle::Particle () : location (), reference_location(), id (0), property_pool(NULL), - properties(PropertyPool::invalid_handle) + properties(PropertyPool::invalid_handle) { } - template - Particle::Particle (const Point &new_location, - const Point &new_reference_location, - const types::particle_index new_id) + + template + Particle::Particle (const Point &location, + const Point &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::invalid_handle) { } - template - Particle::Particle (const Particle &particle) + + + template + Particle::Particle (const Particle &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::invalid_handle) { if (property_pool != NULL) { - const ArrayView my_properties = property_pool->get_properties(properties); + const ArrayView my_properties = property_pool->get_properties(properties); if (my_properties.size() != 0) { - const ArrayView their_properties = particle.get_properties(); + const ArrayView their_properties = particle.get_properties(); std::copy(&their_properties[0],&their_properties[0]+their_properties.size(),&my_properties[0]); } @@ -67,9 +70,10 @@ namespace Particles } - template - Particle::Particle (const void *&data, - PropertyPool &new_property_pool) + + template + Particle::Particle (const void *&data, + PropertyPool &new_property_pool) { const types::particle_index *id_data = static_cast (data); id = *id_data++; @@ -84,16 +88,19 @@ namespace Particles 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 particle_properties = property_pool->get_properties(properties); + const ArrayView particle_properties = property_pool->get_properties(properties); for (unsigned int i = 0; i < particle_properties.size(); ++i) particle_properties[i] = *pdata++; data = static_cast (pdata); } - template - Particle::Particle (Particle &&particle) + + + template + Particle::Particle (Particle &&particle) : location (particle.location), reference_location(particle.reference_location), @@ -101,12 +108,14 @@ namespace Particles property_pool(particle.property_pool), properties(particle.properties) { - particle.properties = PropertyPool::invalid_handle; + particle.properties = PropertyPool::invalid_handle; } - template - Particle & - Particle::operator=(const Particle &particle) + + + template + Particle & + Particle::operator=(const Particle &particle) { if (this != &particle) { @@ -118,23 +127,25 @@ namespace Particles if (property_pool != NULL) { properties = property_pool->allocate_properties_array(); - const ArrayView their_properties = particle.get_properties(); + const ArrayView their_properties = particle.get_properties(); if (their_properties.size() != 0) { - const ArrayView my_properties = property_pool->get_properties(properties); + const ArrayView 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::invalid_handle; } return *this; } - template - Particle & - Particle::operator=(Particle &&particle) + + + template + Particle & + Particle::operator=(Particle &&particle) { if (this != &particle) { @@ -143,21 +154,25 @@ namespace Particles id = particle.id; property_pool = particle.property_pool; properties = particle.properties; - particle.properties = PropertyPool::invalid_handle; + particle.properties = PropertyPool::invalid_handle; } return *this; } - template - Particle::~Particle () + + + template + Particle::~Particle () { - if (properties != PropertyPool::invalid_handle) + if (properties != PropertyPool::invalid_handle) property_pool->deallocate_properties_array(properties); } - template + + + template void - Particle::write_data (void *&data) const + Particle::write_data (void *&data) const { types::particle_index *id_data = static_cast (data); *id_data = id; @@ -172,71 +187,98 @@ namespace Particles for (unsigned int i = 0; i < dim; ++i,++pdata) *pdata = reference_location(i); + // TODO reinterpret cast pdata to PropertyType // Write property data - const ArrayView particle_properties = property_pool->get_properties(properties); + const ArrayView 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 (pdata); } - template + + + template void - Particle::set_location (const Point &new_loc) + Particle::set_location (const Point &new_loc) { location = new_loc; } - template - const Point & - Particle::get_location () const + + + template + const Point & + Particle::get_location () const { return location; } - template + + + template void - Particle::set_reference_location (const Point &new_loc) + Particle::set_reference_location (const Point &new_loc) { reference_location = new_loc; } - template + + + template const Point & - Particle::get_reference_location () const + Particle::get_reference_location () const { return reference_location; } - template + + + template types::particle_index - Particle::get_id () const + Particle::get_id () const { return id; } - template + + + template void - Particle::set_property_pool (PropertyPool &new_property_pool) + Particle::set_property_pool (PropertyPool &new_property_pool) { property_pool = &new_property_pool; } - template + + + template + bool + Particle::has_properties () const + { + return (property_pool != NULL) + && (properties != PropertyPool::invalid_handle); + } + + + + template void - Particle::set_properties (const std::vector &new_properties) + Particle::set_properties (const ArrayView &new_properties) { - if (properties == PropertyPool::invalid_handle) + if (properties == PropertyPool::invalid_handle) properties = property_pool->allocate_properties_array(); - const ArrayView old_properties = property_pool->get_properties(properties); + const ArrayView old_properties = property_pool->get_properties(properties); std::copy(new_properties.begin(),new_properties.end(),&old_properties[0]); } - template - const ArrayView - Particle::get_properties () const + + + template + const ArrayView + Particle::get_properties () const { Assert(property_pool != NULL, ExcInternalError()); @@ -245,9 +287,10 @@ namespace Particles } - template - const ArrayView - Particle::get_properties () + + template + const ArrayView + Particle::get_properties () { Assert(property_pool != NULL, ExcInternalError()); diff --git a/source/particles/particle.inst.in b/source/particles/particle.inst.in index 4b693348c0..9ea1b7927d 100644 --- a/source/particles/particle.inst.in +++ b/source/particles/particle.inst.in @@ -14,11 +14,11 @@ // --------------------------------------------------------------------- -for (deal_II_dimension : DIMENSIONS) +for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS; number : REAL_SCALARS) { namespace Particles \{ template - class Particle ; + class Particle ; \} } diff --git a/source/particles/property_pool.cc b/source/particles/property_pool.cc index 4a85e3fffd..a11d773c86 100644 --- a/source/particles/property_pool.cc +++ b/source/particles/property_pool.cc @@ -14,51 +14,10 @@ // --------------------------------------------------------------------- -#include +#include 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 - PropertyPool::get_properties (const Handle handle) - { - return ArrayView(handle, n_properties); - } - - - - void - PropertyPool::reserve(const std::size_t size) - { - (void)size; - } -} +#include "property_pool.inst" DEAL_II_NAMESPACE_CLOSE diff --git a/source/particles/property_pool.inst.in b/source/particles/property_pool.inst.in new file mode 100644 index 0000000000..974d14254e --- /dev/null +++ b/source/particles/property_pool.inst.in @@ -0,0 +1,24 @@ +// --------------------------------------------------------------------- +// +// 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 ; + \} +} diff --git a/tests/particles/particle_01.cc b/tests/particles/particle_01.cc index 876e03ac16..bfaa707f11 100644 --- a/tests/particles/particle_01.cc +++ b/tests/particles/particle_01.cc @@ -19,8 +19,6 @@ #include "../tests.h" #include -#include -#include template diff --git a/tests/particles/particle_02.cc b/tests/particles/particle_02.cc index 95cc6929f8..e4601726a4 100644 --- a/tests/particles/particle_02.cc +++ b/tests/particles/particle_02.cc @@ -19,8 +19,6 @@ #include "../tests.h" #include -#include -#include template @@ -40,24 +38,24 @@ void test () 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; diff --git a/tests/particles/particle_03.cc b/tests/particles/particle_03.cc index 06d2685f85..9a9e04d2ce 100644 --- a/tests/particles/particle_03.cc +++ b/tests/particles/particle_03.cc @@ -19,8 +19,7 @@ #include "../tests.h" #include -#include -#include +#include template @@ -28,7 +27,7 @@ void test () { { 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; @@ -44,23 +43,28 @@ void test () Particles::Particle<2> particle(position,reference_position,index); particle.set_property_pool(pool); - particle.set_properties(properties); + particle.set_properties(ArrayView(&properties[0],properties.size())); deallog << "Particle properties: " - << std::vector(particle.get_properties().begin(),particle.get_properties().end()) - << std::endl; + << std::vector(particle.get_properties().begin(),particle.get_properties().end()) + << std::endl; const Particles::Particle<2> copy(particle); deallog << "Copy particle properties: " - << std::vector(copy.get_properties().begin(),copy.get_properties().end()) - << std::endl; + << std::vector(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(moved_particle.get_properties().begin(),moved_particle.get_properties().end()) - << std::endl; + << std::vector(moved_particle.get_properties().begin(),moved_particle.get_properties().end()) + << std::endl; + } deallog << "OK" << std::endl; diff --git a/tests/particles/particle_03.output b/tests/particles/particle_03.output index 558d9d51af..1f1af1dc12 100644 --- a/tests/particles/particle_03.output +++ b/tests/particles/particle_03.output @@ -1,5 +1,7 @@ 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 diff --git a/tests/particles/property_pool_01.cc b/tests/particles/property_pool_01.cc new file mode 100644 index 0000000000..e178ca8705 --- /dev/null +++ b/tests/particles/property_pool_01.cc @@ -0,0 +1,51 @@ +// --------------------------------------------------------------------- +// +// 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 +#include +#include + + +template +void test () +{ + { + Particles::PropertyPool pool; + + typename Particles::PropertyPool::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(); + test(); +} diff --git a/tests/particles/property_pool_01.output b/tests/particles/property_pool_01.output new file mode 100644 index 0000000000..fdb9accf31 --- /dev/null +++ b/tests/particles/property_pool_01.output @@ -0,0 +1,5 @@ + +DEAL::Pool properties: 2.50000 +DEAL::OK +DEAL::Pool properties: 2.50000 +DEAL::OK diff --git a/tests/particles/property_pool_02.cc b/tests/particles/property_pool_02.cc new file mode 100644 index 0000000000..f60f69b8e6 --- /dev/null +++ b/tests/particles/property_pool_02.cc @@ -0,0 +1,60 @@ +// --------------------------------------------------------------------- +// +// 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 +#include +#include + + +template +void test () +{ + { + const unsigned int n_properties = 3; + Particles::PropertyPool pool(n_properties); + + typename Particles::PropertyPool::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(); + test(); +} diff --git a/tests/particles/property_pool_02.output b/tests/particles/property_pool_02.output new file mode 100644 index 0000000000..80f767455c --- /dev/null +++ b/tests/particles/property_pool_02.output @@ -0,0 +1,5 @@ + +DEAL::Pool properties: 1.20000 2.50000 2.70000 +DEAL::OK +DEAL::Pool properties: 1.20000 2.50000 2.70000 +DEAL::OK diff --git a/tests/particles/property_pool_03.cc b/tests/particles/property_pool_03.cc new file mode 100644 index 0000000000..48570962f3 --- /dev/null +++ b/tests/particles/property_pool_03.cc @@ -0,0 +1,71 @@ +// --------------------------------------------------------------------- +// +// 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 +#include +#include + + + +struct DataType +{ + double a[2]; + char c; + double *pointer; +}; + +void test () +{ + { + double b[2] = {0.7,1.3}; + + Particles::PropertyPool pool; + + typename Particles::PropertyPool::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(); +} diff --git a/tests/particles/property_pool_03.output b/tests/particles/property_pool_03.output new file mode 100644 index 0000000000..bf8cb587f1 --- /dev/null +++ b/tests/particles/property_pool_03.output @@ -0,0 +1,3 @@ + +DEAL::Pool properties: 2.50000 2.00000 f 1.30000 +DEAL::OK -- 2.39.5