]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use particles with fully distributed tria
authorBruno <blais.bruno@gmail.com>
Wed, 13 Jan 2021 17:48:50 +0000 (12:48 -0500)
committerBruno <blais.bruno@gmail.com>
Wed, 13 Jan 2021 17:48:50 +0000 (12:48 -0500)
source/particles/particle_handler.cc
tests/particles/particle_handler_fully_distributed_01.cc [new file with mode: 0644]
tests/particles/particle_handler_fully_distributed_01.mpirun=2.output [new file with mode: 0644]

index 93bbed5d120978573347b4bb33604770aeddad94..0f5c9b7ed4a24fdd86718cb5890d17c306f7e7db 100644 (file)
@@ -490,7 +490,7 @@ namespace Particles
       AssertDimension(ids.size(), positions.size());
 
     const auto tria =
-      dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim> *>(
+      dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
         &(*triangulation));
     const auto comm =
       (tria != nullptr ? tria->get_communicator() : MPI_COMM_WORLD);
@@ -1799,8 +1799,8 @@ namespace Particles
 
     switch (status)
       {
-        case parallel::distributed::Triangulation<dim, spacedim>::CELL_PERSIST:
-        case parallel::distributed::Triangulation<dim, spacedim>::CELL_REFINE:
+        case parallel::TriangulationBase<dim, spacedim>::CELL_PERSIST:
+        case parallel::TriangulationBase<dim, spacedim>::CELL_REFINE:
           // If the cell persist or is refined store all particles of the
           // current cell.
           {
@@ -1828,7 +1828,7 @@ namespace Particles
           }
           break;
 
-        case parallel::distributed::Triangulation<dim, spacedim>::CELL_COARSEN:
+        case parallel::TriangulationBase<dim, spacedim>::CELL_COARSEN:
           // If this cell is the parent of children that will be coarsened,
           // collect the particles of all children.
           {
@@ -1892,7 +1892,7 @@ namespace Particles
 
     switch (status)
       {
-        case parallel::distributed::Triangulation<dim, spacedim>::CELL_PERSIST:
+        case parallel::TriangulationBase<dim, spacedim>::CELL_PERSIST:
           {
             auto position_hint = particles.end();
             for (const auto &particle : loaded_particles_on_cell)
@@ -1912,7 +1912,7 @@ namespace Particles
           }
           break;
 
-        case parallel::distributed::Triangulation<dim, spacedim>::CELL_COARSEN:
+        case parallel::TriangulationBase<dim, spacedim>::CELL_COARSEN:
           {
             typename std::multimap<internal::LevelInd,
                                    Particle<dim, spacedim>>::iterator
@@ -1938,7 +1938,7 @@ namespace Particles
           }
           break;
 
-        case parallel::distributed::Triangulation<dim, spacedim>::CELL_REFINE:
+        case parallel::TriangulationBase<dim, spacedim>::CELL_REFINE:
           {
             std::vector<
               typename std::multimap<internal::LevelInd,
diff --git a/tests/particles/particle_handler_fully_distributed_01.cc b/tests/particles/particle_handler_fully_distributed_01.cc
new file mode 100644 (file)
index 0000000..2b33298
--- /dev/null
@@ -0,0 +1,147 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 - 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// like particle_handler_04, but for fully distributed triangulations in
+// parallel computations
+
+#include <deal.II/distributed/fully_distributed_tria.h>
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/fe/mapping_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+
+#include <deal.II/particles/particle_handler.h>
+
+#include "../tests.h"
+
+template <int dim, int spacedim>
+void
+test()
+{
+  {
+    const MPI_Comm comm = MPI_COMM_WORLD;
+
+    // create a distributed triangulation
+    parallel::distributed::Triangulation<dim, spacedim> tria_pdt(comm);
+
+    GridGenerator::hyper_cube(tria_pdt);
+    tria_pdt.refine_global(2);
+    MappingQ<dim, spacedim> mapping(1);
+
+
+    // create a fully distributed triangulation
+    parallel::fullydistributed::Triangulation<dim, spacedim> tria_pft(comm);
+
+    // extract relevant information from distributed triangulation
+    auto construction_data = TriangulationDescription::Utilities::
+      create_description_from_triangulation(tria_pdt, comm);
+
+    // actually create triangulation
+    tria_pft.create_triangulation(construction_data);
+
+    //    both processes create a particle handler,
+    //      but only the first creates particles
+    Particles::ParticleHandler<dim, spacedim> particle_handler(tria_pft,
+                                                               mapping);
+
+    if (Utilities::MPI::this_mpi_process(tria_pft.get_communicator()) == 0)
+      {
+        std::vector<Point<spacedim>> position(2);
+        std::vector<Point<dim>>      reference_position(2);
+
+        for (unsigned int i = 0; i < dim; ++i)
+          {
+            position[0](i) = 0.125;
+            position[1](i) = 0.525;
+          }
+
+        Particles::Particle<dim, spacedim> particle1(position[0],
+                                                     reference_position[0],
+                                                     0);
+        Particles::Particle<dim, spacedim> particle2(position[1],
+                                                     reference_position[1],
+                                                     1);
+
+        typename Triangulation<dim, spacedim>::active_cell_iterator cell1(
+          &tria_pft, 2, 0);
+        typename Triangulation<dim, spacedim>::active_cell_iterator cell2(
+          &tria_pft, 2, 0);
+
+        particle_handler.insert_particle(particle1, cell1);
+        particle_handler.insert_particle(particle2, cell2);
+
+        for (const auto &particle : particle_handler)
+          deallog << "Before sort particle id " << particle.get_id()
+                  << " is in cell " << particle.get_surrounding_cell(tria_pft)
+                  << " on process "
+                  << Utilities::MPI::this_mpi_process(
+                       tria_pft.get_communicator())
+                  << std::flush << std::endl;
+      }
+
+
+
+    particle_handler.sort_particles_into_subdomains_and_cells();
+
+    for (const auto &particle : particle_handler)
+      deallog << "After sort particle id " << particle.get_id()
+              << " is in cell " << particle.get_surrounding_cell(tria_pft)
+              << " on process "
+              << Utilities::MPI::this_mpi_process(tria_pft.get_communicator())
+              << std::flush << std::endl;
+
+    // Move all points up by 0.5. This will change cell for particle 1 and will
+    // move particle 2 out of the domain. Note that we need to change the
+    // coordinate dim-1 despite having a spacedim point.
+    Point<spacedim> shift;
+    shift(dim - 1) = 0.5;
+    for (auto &particle : particle_handler)
+      particle.set_location(particle.get_location() + shift);
+
+    particle_handler.sort_particles_into_subdomains_and_cells();
+    for (const auto &particle : particle_handler)
+      deallog << "After shift particle id " << particle.get_id()
+              << " is in cell " << particle.get_surrounding_cell(tria_pft)
+              << " on process "
+              << Utilities::MPI::this_mpi_process(tria_pft.get_communicator())
+              << std::flush << std::endl;
+  }
+
+  deallog << "OK" << std::endl;
+}
+
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  MPILogInitAll all;
+
+  deallog.push("2d/2d");
+  test<2, 2>();
+  deallog.pop();
+  deallog.push("2d/3d");
+  test<2, 3>();
+  deallog.pop();
+  deallog.push("3d/3d");
+  test<3, 3>();
+  deallog.pop();
+}
diff --git a/tests/particles/particle_handler_fully_distributed_01.mpirun=2.output b/tests/particles/particle_handler_fully_distributed_01.mpirun=2.output
new file mode 100644 (file)
index 0000000..01df526
--- /dev/null
@@ -0,0 +1,24 @@
+
+DEAL:0:2d/2d::Before sort particle id 0 is in cell 2.0 on process 0
+DEAL:0:2d/2d::Before sort particle id 1 is in cell 2.0 on process 0
+DEAL:0:2d/2d::After sort particle id 0 is in cell 2.0 on process 0
+DEAL:0:2d/2d::OK
+DEAL:0:2d/3d::Before sort particle id 0 is in cell 2.0 on process 0
+DEAL:0:2d/3d::Before sort particle id 1 is in cell 2.0 on process 0
+DEAL:0:2d/3d::After sort particle id 0 is in cell 2.0 on process 0
+DEAL:0:2d/3d::OK
+DEAL:0:3d/3d::Before sort particle id 0 is in cell 2.0 on process 0
+DEAL:0:3d/3d::Before sort particle id 1 is in cell 2.0 on process 0
+DEAL:0:3d/3d::After sort particle id 0 is in cell 2.0 on process 0
+DEAL:0:3d/3d::OK
+
+DEAL:1:2d/2d::After sort particle id 1 is in cell 2.12 on process 1
+DEAL:1:2d/2d::After shift particle id 0 is in cell 2.8 on process 1
+DEAL:1:2d/2d::OK
+DEAL:1:2d/3d::After sort particle id 1 is in cell 2.12 on process 1
+DEAL:1:2d/3d::After shift particle id 0 is in cell 2.8 on process 1
+DEAL:1:2d/3d::OK
+DEAL:1:3d/3d::After sort particle id 1 is in cell 2.56 on process 1
+DEAL:1:3d/3d::After shift particle id 0 is in cell 2.32 on process 1
+DEAL:1:3d/3d::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.