]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Initial particle implementation
authorRene Gassmoeller <rene.gassmoeller@mailbox.org>
Tue, 8 Aug 2017 20:01:53 +0000 (14:01 -0600)
committerRene Gassmoeller <rene.gassmoeller@mailbox.org>
Fri, 22 Sep 2017 15:59:05 +0000 (09:59 -0600)
include/deal.II/particle/particle.h [new file with mode: 0644]
include/deal.II/particle/property_pool.h [new file with mode: 0644]
source/CMakeLists.txt
source/particle/CMakeLists.txt [new file with mode: 0644]
source/particle/particle.cc [new file with mode: 0644]
source/particle/particle.inst.in [new file with mode: 0644]
source/particle/property_pool.cc [new file with mode: 0644]
tests/particle/CMakeLists.txt [new file with mode: 0644]
tests/particle/particle_01.cc [new file with mode: 0644]
tests/particle/particle_01.output [new file with mode: 0644]

diff --git a/include/deal.II/particle/particle.h b/include/deal.II/particle/particle.h
new file mode 100644 (file)
index 0000000..ab48876
--- /dev/null
@@ -0,0 +1,309 @@
+// ---------------------------------------------------------------------
+//
+// 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__particle_particle_h
+#define dealii__particle_particle_h
+
+#include <deal.II/particle/property_pool.h>
+
+#include <deal.II/base/point.h>
+#include <deal.II/base/types.h>
+#include <deal.II/base/array_view.h>
+
+#include <boost/serialization/vector.hpp>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace Particle
+{
+  /**
+   * A namespace for all type definitions related to particles.
+   */
+  namespace types
+  {
+    /* Type definitions */
+
+#ifdef DEAL_II_WITH_64BIT_INDICES
+    /**
+     * The type used for indices of particles. While in
+     * sequential computations the 4 billion indices of 32-bit unsigned integers
+     * is plenty, parallel computations using hundreds of processes can overflow
+     * this number and we need a bigger index space. We here utilize the same
+     * build variable that controls the dof indices because the number
+     * of degrees of freedom and the number of particles are typically on the same
+     * order of magnitude.
+     *
+     * The data type always indicates an unsigned integer type.
+     */
+    typedef unsigned long long int particle_index;
+
+    /**
+     * An identifier that denotes the MPI type associated with
+     * types::global_dof_index.
+     */
+#  define ASPECT_PARTICLE_INDEX_MPI_TYPE MPI_UNSIGNED_LONG_LONG
+#else
+    /**
+     * The type used for indices of particles. While in
+     * sequential computations the 4 billion indices of 32-bit unsigned integers
+     * is plenty, parallel computations using hundreds of processes can overflow
+     * this number and we need a bigger index space. We here utilize the same
+     * build variable that controls the dof indices because the number
+     * of degrees of freedom and the number of particles are typically on the same
+     * order of magnitude.
+     *
+     * The data type always indicates an unsigned integer type.
+     */
+    typedef unsigned int particle_index;
+
+    /**
+     * An identifier that denotes the MPI type associated with
+     * types::global_dof_index.
+     */
+#  define ASPECT_PARTICLE_INDEX_MPI_TYPE MPI_UNSIGNED
+#endif
+  }
+
+  /**
+   * 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.
+   *
+   *
+   * @ingroup Particle
+   *
+   */
+  template <int dim>
+  class Particle
+  {
+  public:
+    /**
+     * Empty constructor for Particle, creates a particle at the
+     * origin.
+     */
+    Particle ();
+
+    /**
+     * Constructor for Particle, creates a particle with the specified
+     * ID at the specified location. Note that Aspect
+     * does not check for duplicate particle IDs so the generator must
+     * make sure the IDs are unique over all processes.
+     *
+     * @param[in] new_location Initial location of particle.
+     * @param[in] new_reference_location Initial location of the particle
+     * in the coordinate system of the reference cell.
+     * @param[in] new_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);
+
+    /**
+     * Copy-Constructor for Particle, creates a particle with exactly the
+     * state of the input argument. Note that since each particle has a
+     * handle for a certain piece of the property memory, and is responsible
+     * 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);
+
+    /**
+     * Constructor for Particle, creates a particle from a data vector.
+     * This constructor is usually called after sending a particle to a
+     * different process.
+     *
+     * @param[in,out] begin_data A pointer to a memory location from which
+     * to read the information that completely describes a particle. This
+     * class then de-serializes its data from this memory location and
+     * advance the pointer accordingly.
+     *
+     * @param[in,out] new_property_pool A property pool that is used to
+     * allocate the property data used by this particle.
+     */
+    Particle (const void *&begin_data,
+              PropertyPool &new_property_pool);
+
+    /**
+     * Move constructor for Particle, creates a particle from an existing
+     * one by stealing its state.
+     */
+    Particle (Particle<dim> &&particle);
+
+    /**
+     * Copy assignment operator.
+     */
+    Particle<dim> &operator=(const Particle<dim> &particle);
+
+    /**
+     * Move assignment operator.
+     */
+    Particle<dim> &operator=(Particle<dim> &&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).
+     */
+    ~Particle ();
+
+    /**
+     * 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, const unsigned int data_size);
+     *
+     * @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_loc The new location for this particle.
+     */
+    void
+    set_location (const Point<dim> &new_loc);
+
+    /**
+     * Get the location of this particle.
+     *
+     * @return The location of this particle.
+     */
+    const Point<dim> &
+    get_location () const;
+
+    /**
+     * Set the reference 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 reference location for this particle.
+     */
+    void
+    set_reference_location (const Point<dim> &new_loc);
+
+    /**
+     * Get the reference location of this particle in its current cell.
+     *
+     * @return The reference location of this particle.
+     */
+    const Point<dim> &
+    get_reference_location () const;
+
+    /**
+     * Get the ID number of this particle.
+     *
+     * @return The id of this particle.
+     */
+    types::particle_index
+    get_id () const;
+
+    /**
+     * Tell the particle where to store its properties (even if it does not
+     * own properties). Usually this is only done once per particle, but
+     * since the particle generator does not know about the properties
+     * we want to do it not at construction time. Another use for this
+     * function is after particle transfer to a new process.
+     */
+    void
+    set_property_pool(PropertyPool &property_pool);
+
+    /**
+     * 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;
+
+    /**
+     * Serialize the contents of this class.
+     */
+    template <class Archive>
+    void serialize (Archive &ar, const unsigned int version);
+
+  private:
+    /**
+     * Current particle location.
+     */
+    Point<dim>             location;
+
+    /**
+     * Current particle location in the reference cell.
+     * Storing this reduces the number of times we need to compute this
+     * location, which takes a significant amount of computing time.
+     */
+    Point<dim>             reference_location;
+
+    /**
+     * Globally unique ID of particle.
+     */
+    types::particle_index  id;
+
+    /**
+     * A pointer to the property pool. Necessary to translate from the
+     * handle to the actual memory locations.
+     */
+    PropertyPool *property_pool;
+
+    /**
+     * A handle to all particle properties
+     */
+    PropertyPool::Handle properties;
+  };
+
+  /* -------------------------- inline and template functions ---------------------- */
+
+  template <int dim>
+  template <class Archive>
+  void Particle<dim>::serialize (Archive &ar, const unsigned int)
+  {
+    ar &location
+    & id
+    & properties
+    ;
+  }
+}
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
+
diff --git a/include/deal.II/particle/property_pool.h b/include/deal.II/particle/property_pool.h
new file mode 100644 (file)
index 0000000..99dc52c
--- /dev/null
@@ -0,0 +1,88 @@
+// ---------------------------------------------------------------------
+//
+// 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__particle_property_pool_h
+#define dealii__particle_property_pool_h
+
+#include <deal.II/base/array_view.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace Particle
+{
+  /**
+   * 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.
+   */
+  class PropertyPool
+  {
+  public:
+    /**
+     * Typedef for the handle that is returned to the particles, and that
+     * uniquely identifies the slot of memory that is reserved for this
+     * particle.
+     */
+    typedef double *Handle;
+
+    /**
+     * Define a default (invalid) value for handles.
+     */
+    static const Handle invalid_handle;
+
+    /**
+     * Constructor. Stores the number of properties per reserved slot.
+     */
+    PropertyPool (const unsigned int n_properties_per_slot);
+
+    /**
+     * Returns a new handle that allows accessing the reserved block
+     * of memory.
+     */
+    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.
+     */
+    void deallocate_properties_array (const Handle handle);
+
+    /**
+     * Return an ArrayView to the properties that correspond to the given
+     * handle @p handle.
+     */
+    ArrayView<double> get_properties (const Handle handle);
+
+    /**
+     * Reserves the dynamic memory needed for storing the properties of
+     * @p size particles.
+     */
+    void reserve(const std::size_t size);
+
+  private:
+    /**
+     * The number of properties that are reserved per particle.
+     */
+    const unsigned int n_properties;
+  };
+
+}
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index 6b772be0141b287ffa5956ba46a70ee2d97267f6..e21a35381102f161dbb62ab671e5958193b73df2 100644 (file)
@@ -48,6 +48,7 @@ ADD_SUBDIRECTORY(integrators)
 ADD_SUBDIRECTORY(matrix_free)
 ADD_SUBDIRECTORY(meshworker)
 ADD_SUBDIRECTORY(opencascade)
+ADD_SUBDIRECTORY(particle)
 ADD_SUBDIRECTORY(physics)
 ADD_SUBDIRECTORY(non_matching)
 ADD_SUBDIRECTORY(sundials)
diff --git a/source/particle/CMakeLists.txt b/source/particle/CMakeLists.txt
new file mode 100644 (file)
index 0000000..99db585
--- /dev/null
@@ -0,0 +1,33 @@
+## ---------------------------------------------------------------------
+##
+## Copyright (C) 2016 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_DIRECTORIES(BEFORE ${CMAKE_CURRENT_BINARY_DIR})
+
+SET(_src
+  particle.cc
+  property_pool.cc
+  )
+
+SET(_inst
+  particle.inst.in
+  )
+
+FILE(GLOB _header
+  ${CMAKE_SOURCE_DIR}/include/deal.II/particle/*.h
+  )
+
+DEAL_II_ADD_LIBRARY(obj_particle OBJECT ${_src} ${_header} ${_inst})
+EXPAND_INSTANTIATIONS(obj_particle "${_inst}")
diff --git a/source/particle/particle.cc b/source/particle/particle.cc
new file mode 100644 (file)
index 0000000..f842b41
--- /dev/null
@@ -0,0 +1,268 @@
+// ---------------------------------------------------------------------
+//
+// 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/particle/particle.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace Particle
+{
+  template <int dim>
+  Particle<dim>::Particle ()
+    :
+    location (),
+    reference_location(),
+    id (0),
+    property_pool(NULL),
+    properties(PropertyPool::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)
+    :
+    location (new_location),
+    reference_location (new_reference_location),
+    id (new_id),
+    property_pool(NULL),
+    properties (PropertyPool::invalid_handle)
+  {
+  }
+
+  template <int dim>
+  Particle<dim>::Particle (const Particle<dim> &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)
+  {
+    if (property_pool != NULL)
+      {
+        const ArrayView<double> my_properties = property_pool->get_properties(properties);
+
+        if (my_properties.size() != 0)
+          {
+            const ArrayView<const double> 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)
+  {
+    const types::particle_index *id_data = static_cast<const types::particle_index *> (data);
+    id = *id_data++;
+    const double *pdata = reinterpret_cast<const double *> (id_data);
+
+    for (unsigned int i = 0; i < dim; ++i)
+      location(i) = *pdata++;
+
+    for (unsigned int i = 0; i < dim; ++i)
+      reference_location(i) = *pdata++;
+
+    property_pool = &new_property_pool;
+    properties = property_pool->allocate_properties_array();
+
+    // See if there are properties to load
+    const ArrayView<double> 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);
+  }
+
+#ifdef DEAL_II_WITH_CXX11
+
+  template <int dim>
+  Particle<dim>::Particle (Particle<dim> &&particle)
+    :
+    location (particle.location),
+    reference_location(particle.reference_location),
+    id (particle.id),
+    property_pool(particle.property_pool),
+    properties(particle.properties)
+  {
+    particle.properties = PropertyPool::invalid_handle;
+  }
+
+  template <int dim>
+  Particle<dim> &
+  Particle<dim>::operator=(const Particle<dim> &particle)
+  {
+    if (this != &particle)
+      {
+        location = particle.location;
+        reference_location = particle.reference_location;
+        id = particle.id;
+        property_pool = particle.property_pool;
+
+        if (property_pool != NULL)
+          {
+            properties = property_pool->allocate_properties_array();
+            const ArrayView<const double> their_properties = particle.get_properties();
+
+            if (their_properties.size() != 0)
+              {
+                const ArrayView<double> 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;
+      }
+    return *this;
+  }
+
+  template <int dim>
+  Particle<dim> &
+  Particle<dim>::operator=(Particle<dim> &&particle)
+  {
+    if (this != &particle)
+      {
+        location = particle.location;
+        reference_location = particle.reference_location;
+        id = particle.id;
+        property_pool = particle.property_pool;
+        properties = particle.properties;
+        particle.properties = PropertyPool::invalid_handle;
+      }
+    return *this;
+  }
+#endif
+
+  template <int dim>
+  Particle<dim>::~Particle ()
+  {
+    if (properties != PropertyPool::invalid_handle)
+      property_pool->deallocate_properties_array(properties);
+  }
+
+  template <int dim>
+  void
+  Particle<dim>::write_data (void *&data) const
+  {
+    types::particle_index *id_data  = static_cast<types::particle_index *> (data);
+    *id_data = id;
+    ++id_data;
+    double *pdata = reinterpret_cast<double *> (id_data);
+
+    // Write location data
+    for (unsigned int i = 0; i < dim; ++i,++pdata)
+      *pdata = location(i);
+
+    // Write reference location data
+    for (unsigned int i = 0; i < dim; ++i,++pdata)
+      *pdata = reference_location(i);
+
+    // Write property data
+    const ArrayView<double> 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>
+  void
+  Particle<dim>::set_location (const Point<dim> &new_loc)
+  {
+    location = new_loc;
+  }
+
+  template <int dim>
+  const Point<dim> &
+  Particle<dim>::get_location () const
+  {
+    return location;
+  }
+
+  template <int dim>
+  void
+  Particle<dim>::set_reference_location (const Point<dim> &new_loc)
+  {
+    reference_location = new_loc;
+  }
+
+  template <int dim>
+  const Point<dim> &
+  Particle<dim>::get_reference_location () const
+  {
+    return reference_location;
+  }
+
+  template <int dim>
+  types::particle_index
+  Particle<dim>::get_id () const
+  {
+    return id;
+  }
+
+  template <int dim>
+  void
+  Particle<dim>::set_property_pool (PropertyPool &new_property_pool)
+  {
+    property_pool = &new_property_pool;
+  }
+
+  template <int dim>
+  void
+  Particle<dim>::set_properties (const std::vector<double> &new_properties)
+  {
+    if (properties == PropertyPool::invalid_handle)
+      properties = property_pool->allocate_properties_array();
+
+    const ArrayView<double> 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
+  {
+    Assert(property_pool != NULL,
+           ExcInternalError());
+
+    return property_pool->get_properties(properties);
+  }
+
+
+  template <int dim>
+  const ArrayView<double>
+  Particle<dim>::get_properties ()
+  {
+    Assert(property_pool != NULL,
+           ExcInternalError());
+
+    return property_pool->get_properties(properties);
+  }
+}
+
+DEAL_II_NAMESPACE_CLOSE
+
+DEAL_II_NAMESPACE_OPEN
+
+#include "particle.inst"
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/source/particle/particle.inst.in b/source/particle/particle.inst.in
new file mode 100644 (file)
index 0000000..0c888d7
--- /dev/null
@@ -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 (deal_II_dimension : DIMENSIONS)
+{
+    namespace Particle
+    \{
+    template
+    class Particle <deal_II_dimension>;
+    \}
+}
diff --git a/source/particle/property_pool.cc b/source/particle/property_pool.cc
new file mode 100644 (file)
index 0000000..1aa1ada
--- /dev/null
@@ -0,0 +1,64 @@
+// ---------------------------------------------------------------------
+//
+// 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/particle/property_pool.h>
+#include <deal.II/particle/particle.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace Particle
+{
+  const 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;
+  }
+}
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/tests/particle/CMakeLists.txt b/tests/particle/CMakeLists.txt
new file mode 100644 (file)
index 0000000..d33eafa
--- /dev/null
@@ -0,0 +1,4 @@
+CMAKE_MINIMUM_REQUIRED(VERSION 2.8.12)
+INCLUDE(../setup_testsubproject.cmake)
+PROJECT(testsuite CXX)
+DEAL_II_PICKUP_TESTS()
diff --git a/tests/particle/particle_01.cc b/tests/particle/particle_01.cc
new file mode 100644 (file)
index 0000000..fa8a009
--- /dev/null
@@ -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 and destruction of particles
+
+#include "../tests.h"
+#include <deal.II/particle/particle.h>
+#include <fstream>
+#include <iomanip>
+
+
+template <int dim>
+void test ()
+{
+  {
+    Particle::Particle<dim> particle;
+
+    deallog << "Particle location: " << particle.get_location() << std::endl;
+
+    Point<dim> position;
+    position(0) = 1.0;
+    particle.set_location(position);
+
+    deallog << "Particle location: " << particle.get_location() << std::endl;
+  }
+
+  deallog << "OK" << std::endl;
+}
+
+
+
+int main ()
+{
+  initlog();
+  test<2>();
+  test<3>();
+}
diff --git a/tests/particle/particle_01.output b/tests/particle/particle_01.output
new file mode 100644 (file)
index 0000000..cbfa261
--- /dev/null
@@ -0,0 +1,7 @@
+
+DEAL::Particle location: 0.00000 0.00000
+DEAL::Particle location: 1.00000 0.00000
+DEAL::OK
+DEAL::Particle location: 0.00000 0.00000 0.00000
+DEAL::Particle location: 1.00000 0.00000 0.00000
+DEAL::OK

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.