]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add test
authorRene Gassmoeller <rene.gassmoeller@mailbox.org>
Thu, 2 Nov 2017 23:30:24 +0000 (17:30 -0600)
committerRene Gassmoeller <rene.gassmoeller@mailbox.org>
Thu, 9 Nov 2017 16:47:23 +0000 (09:47 -0700)
include/deal.II/particles/particle_handler.h
source/particles/particle_handler.cc
tests/particles/particle_handler_01.cc
tests/particles/particle_handler_01.output
tests/particles/particle_handler_02.cc [new file with mode: 0644]
tests/particles/particle_handler_02.output [new file with mode: 0644]

index 271fd2fea495d39d55005751df234c84f5434f25..15df27dce6bbe4bc77bf232554e6af597d35663d 100644 (file)
@@ -211,6 +211,14 @@ namespace Particles
      */
     types::particle_index n_global_particles() const;
 
+    /**
+     * Return the maximum number of particles per cell the last
+     * time the update_n_global_particles() function was called.
+     *
+     * @return Maximum number of particles in one cell in simulation.
+     */
+    types::particle_index n_global_max_particles_per_cell() const;
+
     /**
      * Return the number of particles in the local part of the
      * triangulation.
index e57a824e8ff2e4ec86458ed4ec5debb5008671b1..53ce6ee87bff308cff62bb87ee49db810139eae8 100644 (file)
@@ -111,14 +111,17 @@ namespace Particles
     types::particle_index locally_highest_index = 0;
     unsigned int local_max_particles_per_cell = 0;
     unsigned int current_particles_per_cell = 0;
-    typename Triangulation<dim,spacedim>::active_cell_iterator current_cell;
+    typename Triangulation<dim,spacedim>::active_cell_iterator current_cell = triangulation->begin_active();
 
     for (particle_iterator particle = begin(); particle != end(); ++particle)
       {
         locally_highest_index = std::max(locally_highest_index,particle->get_id());
 
         if (particle->get_surrounding_cell(*triangulation) != current_cell)
-          current_particles_per_cell = 0;
+          {
+            current_particles_per_cell = 0;
+            current_cell = particle->get_surrounding_cell(*triangulation);
+          }
 
         ++current_particles_per_cell;
         local_max_particles_per_cell = std::max(local_max_particles_per_cell,
@@ -246,6 +249,15 @@ namespace Particles
 
 
 
+  template <int dim,int spacedim>
+  types::particle_index
+  ParticleHandler<dim,spacedim>::n_global_max_particles_per_cell() const
+  {
+    return global_max_particles_per_cell;
+  }
+
+
+
   template <int dim,int spacedim>
   types::particle_index
   ParticleHandler<dim,spacedim>::n_locally_owned_particles() const
index 24557329e5f45678fd40a75b4d7a5098871ebbac..0bb3d527e4e627d8e85e01e4fd16e2223c921a27 100644 (file)
@@ -36,7 +36,6 @@ void test ()
 
     Particles::ParticleHandler<dim,spacedim> particle_handler(tr,mapping);
 
-    Particles::Particle<dim,spacedim> particle;
 
     Point<spacedim> position;
     position(0) = 0.3;
@@ -50,10 +49,9 @@ void test ()
     if (dim>1)
       reference_position(1) = 0.4;
     if (dim>2)
-      reference_position(1) = 0.6;
+      reference_position(2) = 0.6;
 
-    particle.set_location(position);
-    particle.set_reference_location(reference_position);
+    Particles::Particle<dim,spacedim> particle(position,reference_position,7);
     deallog << "Particle location: " << particle.get_location() << std::endl;
 
 
@@ -61,7 +59,7 @@ void test ()
       GridTools::find_active_cell_around_point(mapping, tr, particle.get_location());
 
     particle_handler.insert_particle(particle,cell_position.first);
-    particle_handler.update_cache();
+    particle_handler.update_cached_numbers();
 
     deallog << "Particle number: " << particle_handler.n_global_particles() << std::endl;
 
index 2cf7a2752210d080840b4412613b032cc5b15b01..044ef21d4e73ab4b6a0eb3292d3909cec4a02f15 100644 (file)
@@ -12,5 +12,5 @@ DEAL::OK
 DEAL::Particle location: 0.300000 0.500000 0.700000
 DEAL::Particle number: 1
 DEAL::Particle location: 0.300000 0.500000 0.700000
-DEAL::Particle reference location: 0.200000 0.600000 0.00000
+DEAL::Particle reference location: 0.200000 0.400000 0.600000
 DEAL::OK
diff --git a/tests/particles/particle_handler_02.cc b/tests/particles/particle_handler_02.cc
new file mode 100644 (file)
index 0000000..6d94ddb
--- /dev/null
@@ -0,0 +1,99 @@
+// ---------------------------------------------------------------------
+//
+// 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 particle within the particle handler class.
+
+#include "../tests.h"
+#include <deal.II/particles/particle_handler.h>
+
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/fe/mapping_q.h>
+
+template <int dim, int spacedim>
+void test ()
+{
+  {
+    parallel::distributed::Triangulation<dim,spacedim> tr(MPI_COMM_WORLD);
+
+    GridGenerator::hyper_cube(tr);
+    tr.refine_global(1);
+    MappingQ<dim,spacedim> mapping(1);
+
+    Particles::ParticleHandler<dim,spacedim> particle_handler(tr,mapping);
+
+
+    Point<spacedim> position;
+    position(0) = 0.3;
+    if (spacedim>1)
+      position(1) = 0.5;
+    if (spacedim>2)
+      position(2) = 0.7;
+
+    Point<dim> reference_position;
+    reference_position(0) = 0.2;
+    if (dim>1)
+      reference_position(1) = 0.4;
+    if (dim>2)
+      reference_position(2) = 0.6;
+
+    Particles::Particle<dim,spacedim> particle(position,reference_position,7);
+    deallog << "Particle location: " << particle.get_location() << std::endl;
+
+
+    std::pair<typename parallel::distributed::Triangulation<dim, spacedim>::active_cell_iterator, Point<dim> > cell_position =
+      GridTools::find_active_cell_around_point(mapping, tr, particle.get_location());
+
+    particle_handler.insert_particle(particle,cell_position.first);
+    particle_handler.insert_particle(particle,cell_position.first);
+
+    position(0) = 0.7;
+    Particles::Particle<dim,spacedim> particle2(position,reference_position,9);
+
+    cell_position = GridTools::find_active_cell_around_point(mapping, tr, particle2.get_location());
+    particle_handler.insert_particle(particle2,cell_position.first);
+
+    particle_handler.update_cached_numbers();
+
+    deallog << "Particle number: " << particle_handler.n_global_particles() << std::endl;
+    deallog << "Next free particle index: " << particle_handler.get_next_free_particle_index() << std::endl;
+    deallog << "Max particles per cell: " << particle_handler.n_global_max_particles_per_cell() << std::endl;
+
+    for (auto particle = particle_handler.begin(); particle != particle_handler.end(); ++particle)
+      {
+        deallog << "Particle location: " << particle->get_location() << std::endl;
+        deallog << "Particle reference location: " << particle->get_reference_location() << std::endl;
+      }
+  }
+
+  deallog << "OK" << std::endl;
+}
+
+
+
+int main (int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);
+
+  initlog();
+
+  test<2,2>();
+  test<2,3>();
+
+  test<3,3>();
+}
diff --git a/tests/particles/particle_handler_02.output b/tests/particles/particle_handler_02.output
new file mode 100644 (file)
index 0000000..4238f03
--- /dev/null
@@ -0,0 +1,34 @@
+
+DEAL::Particle location: 0.300000 0.500000
+DEAL::Particle number: 3
+DEAL::Next free particle index: 10
+DEAL::Max particles per cell: 2
+DEAL::Particle location: 0.300000 0.500000
+DEAL::Particle reference location: 0.200000 0.400000
+DEAL::Particle location: 0.300000 0.500000
+DEAL::Particle reference location: 0.200000 0.400000
+DEAL::Particle location: 0.700000 0.500000
+DEAL::Particle reference location: 0.200000 0.400000
+DEAL::OK
+DEAL::Particle location: 0.300000 0.500000 0.700000
+DEAL::Particle number: 3
+DEAL::Next free particle index: 10
+DEAL::Max particles per cell: 2
+DEAL::Particle location: 0.300000 0.500000 0.700000
+DEAL::Particle reference location: 0.200000 0.400000
+DEAL::Particle location: 0.300000 0.500000 0.700000
+DEAL::Particle reference location: 0.200000 0.400000
+DEAL::Particle location: 0.700000 0.500000 0.700000
+DEAL::Particle reference location: 0.200000 0.400000
+DEAL::OK
+DEAL::Particle location: 0.300000 0.500000 0.700000
+DEAL::Particle number: 3
+DEAL::Next free particle index: 10
+DEAL::Max particles per cell: 2
+DEAL::Particle location: 0.300000 0.500000 0.700000
+DEAL::Particle reference location: 0.200000 0.400000 0.600000
+DEAL::Particle location: 0.300000 0.500000 0.700000
+DEAL::Particle reference location: 0.200000 0.400000 0.600000
+DEAL::Particle location: 0.700000 0.500000 0.700000
+DEAL::Particle reference location: 0.200000 0.400000 0.600000
+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.